{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Clustering using GPU Accelerated DBSCAN in RAPIDS\n",
    "#### By Paul Hendricks\n",
    "-------\n",
    "\n",
    "While the world’s data doubles each year, CPU computing has hit a brick wall with the end of Moore’s law. For the same reasons, scientific computing and deep learning has turned to NVIDIA GPU acceleration, data analytics and machine learning where GPU acceleration is ideal. \n",
    "\n",
    "NVIDIA created RAPIDS – an open-source data analytics and machine learning acceleration platform that leverages GPUs to accelerate computations. RAPIDS is based on Python, has pandas-like and Scikit-Learn-like interfaces, is built on Apache Arrow in-memory data format, and can scale from 1 to multi-GPU to multi-nodes. RAPIDS integrates easily into the world’s most popular data science Python-based workflows. RAPIDS accelerates data science end-to-end – from data prep, to machine learning, to deep learning. And through Arrow, Spark users can easily move data into the RAPIDS platform for acceleration.\n",
    "\n",
    "In this notebook, we will also show how to use DBSCAN - a popular clustering algorithm - and how to use the GPU accelerated implementation of this algorithm in RAPIDS.\n",
    "\n",
    "**Table of Contents**\n",
    "\n",
    "* Clustering with DBSCAN\n",
    "* Setup\n",
    "* Generating Data\n",
    "* K Means and Agglomerative Clustering\n",
    "* Clustering using DBSCAN\n",
    "* Accelerating DBSCAN with RAPIDS\n",
    "* Benchmarking: Comparing GPU and CPU\n",
    "* Conclusion\n",
    "\n",
    "Before going any further, let's make sure we have access to `matplotlib`, a popular Python library for data visualization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "\n",
    "try:\n",
    "    import matplotlib\n",
    "except ModuleNotFoundError:\n",
    "    os.system('conda install -y matplotlib')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clustering with DBSCAN\n",
    "\n",
    "Clustering is an important technique for helping data scientists partition data, especially when that data doesn't have labels or annotations associated with it. Since these data often don't have labels, clustering is often described as an unsupervised learning technique. While there are many different algorithms that partition data into unique clusters, we will show in this notebook how in certain cases the DBSCAN algorithm can do a better job of clustering than traditional algorithms such as K Means or Agglomerative Clustering. \n",
    "\n",
    "We will also show how to cluster data with DBSCAN in NVIDIA RAPIDS – an open-source data analytics and machine learning acceleration platform that leverages GPUs to accelerate computations. RAPIDS is based on Python, has pandas-like and Scikit-Learn-like interfaces, is built on Apache Arrow in-memory data format, and can scale from 1 to multi-GPU to multi-nodes. We will see that porting this example from CPU to GPU is trivial and that we can experience massive performance gains by doing so."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setup\n",
    "\n",
    "This notebook was tested using the `nvcr.io/nvidia/rapidsai/rapidsai:0.5-cuda10.0-runtime-ubuntu18.04-gcc7-py3.7` Docker container from [NVIDIA GPU Cloud](https://ngc.nvidia.com) and run on the NVIDIA Tesla V100 GPU. Please be aware that your system may be different and you may need to modify the code or install packages to run the below examples. \n",
    "\n",
    "If you think you have found a bug or an error, please file an issue here: https://github.com/rapidsai/notebooks/issues\n",
    "\n",
    "Before we begin, let's check out our hardware setup by running the `nvidia-smi` command."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Tue May  7 00:34:27 2019       \n",
      "+-----------------------------------------------------------------------------+\n",
      "| NVIDIA-SMI 418.39       Driver Version: 418.39       CUDA Version: 10.1     |\n",
      "|-------------------------------+----------------------+----------------------+\n",
      "| GPU  Name        Persistence-M| Bus-Id        Disp.A | Volatile Uncorr. ECC |\n",
      "| Fan  Temp  Perf  Pwr:Usage/Cap|         Memory-Usage | GPU-Util  Compute M. |\n",
      "|===============================+======================+======================|\n",
      "|   0  Quadro GV100        Off  | 00000000:15:00.0 Off |                  Off |\n",
      "| 29%   40C    P2    26W / 250W |  10149MiB / 32478MiB |      0%      Default |\n",
      "+-------------------------------+----------------------+----------------------+\n",
      "|   1  Quadro GV100        Off  | 00000000:2D:00.0  On |                  Off |\n",
      "| 33%   46C    P0    29W / 250W |    260MiB / 32470MiB |     24%      Default |\n",
      "+-------------------------------+----------------------+----------------------+\n",
      "                                                                               \n",
      "+-----------------------------------------------------------------------------+\n",
      "| Processes:                                                       GPU Memory |\n",
      "|  GPU       PID   Type   Process name                             Usage      |\n",
      "|=============================================================================|\n",
      "+-----------------------------------------------------------------------------+\n"
     ]
    }
   ],
   "source": [
    "!nvidia-smi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's see what CUDA version we have:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "nvcc: NVIDIA (R) Cuda compiler driver\n",
      "Copyright (c) 2005-2018 NVIDIA Corporation\n",
      "Built on Sat_Aug_25_21:08:01_CDT_2018\n",
      "Cuda compilation tools, release 10.0, V10.0.130\n"
     ]
    }
   ],
   "source": [
    "!nvcc --version"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's load some helper functions from `matplotlib` and configure the Jupyter Notebook for visualization."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from matplotlib.colors import ListedColormap\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generating Data\n",
    "\n",
    "We'll generate some fake data using the `make_moons` function from the `sklearn.datasets` module. This function generates data points from two equations, each describing a half circle with a unique center. Since each data point is generated by one of these two equations, the cluster each data point belongs to is clear. The ideal clustering algorithm will identify two clusters and associate each data point with the equation that generated it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(100, 2)\n"
     ]
    }
   ],
   "source": [
    "from sklearn.datasets import make_moons\n",
    "\n",
    "X, y = make_moons(n_samples=int(1e2), noise=0.05, random_state=0)\n",
    "print(X.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's visualize our data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3X+QXWd93/H3R+ulXmPwmliAtZZidUaVMXWM4MZ2Rm1qDI6MSZBCcGsg4ccko3EnzpAZRkUkMwlMMmVTTdPQ8sOjMS4wOHgc2xECC1SwTGlJTb1CFkaWBaqd2Fq5ePkhJ7HV0Ur69o+9K1/dPffnOffe8+PzmtHsvec+e8+5e3TP9zzP832eRxGBmZlZ3iwb9QGYmZklcYAyM7NccoAyM7NccoAyM7NccoAyM7NccoAyM7NccoAyM7NccoAyM7NccoAyM7NcOmfUB9DORRddFJdeeumoD8PMzDK0d+/eH0fE8k7lch2gLr30UmZmZkZ9GGZmliFJf9dNOTfxmZlZLjlAmZlZLjlAmZlZLjlAmZlZLjlAmZlZLjlAmZlZLjlAmZlZLjlAmZlZLuV6oK5ZGjv2zbJt9yGOHjvOiskJtmxYy6Z1U6M+LDPrkgOU5U4WgWXHvlk+fN+jHJ8/BcDsseN8+L5HARykzArCTXyWK4uBZfbYcYIXA8uOfbM9vc+23YfOBKdFx+dPsW33oQyP1swGyTUoy5V2gaWXms/RY8fbbnfzn1n+uQZludIpsHRrxeREy+1Z1dLMbLAcoCxX2gWWXmzZsJaJ8bGztk2Mj7Flw1o3/5kVhAOU5Uq7wNKLTeum+Njbr2BqcgIBU5MTfOztV7Bp3VRmtTQzG6xM+qAk3QH8KvBsRPzzhNcFfBy4EXgBeF9EfDeLfVu5LPYDZdE/tGndVOLvrZicYDYhGPVaS8uK+8PMkmWVJPFZ4BPA51u8/hZgTf3f1cCn6z+t5Lq9+A7zIr1lw9qzUtChv1paFpwOb9ZaJgEqIr4l6dI2RTYCn4+IAB6SNCnp4oh4Jov9Wz51e/Ed9kW631raIIJoVlmLZmU0rDTzKeDphudH6tuWBChJm4HNAKtWrRrKwdlgdHvxHcVFulXzXyuDCqLuDzNrbVgBSgnbIqlgRGwHtgPUarXEMjY4WdYSur34ZnWRHmQz4aCCaN76w8zyZFhZfEeAlQ3PLwGODmnf1qWsxwd1mzKeRWr5oMc2Daqmk1XWolkZDStA7QTeowXXAM+5/yl/sh4f1O3FN4uL9KDHNmU1PqtZu3R4s6rLKs38i8C1wEWSjgB/DIwDRMRtwC4WUswPs5Bm/v4s9mv9S2oOy7qW0G0yQhap5YPuyxlk5l+v/WFmVaGFxLp8qtVqMTMzM+rDKJ3mDn9YuNj+k3OWcez4/JLyU5MTfHvrdcM8xJ6tn96T2JfT67G368fyeCWzbEjaGxG1TuU8WWwFtWoOO3d8GRPjY7kYH9SrLGo4nTL1XNMxGy4HqApq1ex17IV5/tO/eV0hawlZNBP2kqnXXJt642XLefDxucL93czyzAGqgibPG+dnLyxtyps8b7zQtYS0x95tP1ZSTesLDz115nXPBmGWDU8WW0Gtuh1z3B05FN1m6iXVtJoNa3b0HftmWT+9h9Vb72f99B4vGWKl4gBVQc8lJEK0214V3aa7d5sZOOjZILyulZWdA1QFXTAxnri96rMXdDsmqdu/06D/nl7XysrOfVAVs2PfLM+fOLlk+/gyFSJbb9C66cdKyhhsNozsR8/jZ2XnGlTFbNt9iPlTSzubzj/3HHfodymppvWb16wa+mwQg5rdwiwvXIOqmHYp5ta9PGQ75mldK7NBcA2qYnzXXR6ex8/KzjWoihn1XbenC8pWHmpyZoPiAFUxWcy40C8vb25mvXCAqqBR3XV7eXMz64UDVEnlsSnNadFm1gsnSZRQXmcYcIKGmfXCAaqE8jrDgJc3N7NeuImvhPLalDbKBA07Wx6bgM2aOUCV0IrJicTVZfPQlOa06NFzNqUVRSZNfJJukHRI0mFJWxNev0DSlyXtl3RA0vuz2K8lc1OatZPXJmCzZqlrUJLGgE8C1wNHgIcl7YyIxxqK/S7wWET8mqTlwCFJd0bEibT7r5pummbclGbt5LUJ2KxZFk18VwGHI+IJAEl3ARuBxgAVwMskCTgf+CmwdEpta6uXphk3pVkreW4CNmuURRPfFPB0w/Mj9W2NPgG8BjgKPAp8ICJOZ7DvSnHTjGXBTcBWFFkEKCVsa17PYQPwCLACeB3wCUkvT3wzabOkGUkzc3NzGRxeebhpxrLgSWatKLJo4jsCrGx4fgkLNaVG7wemIyKAw5KeBC4D/nfzm0XEdmA7QK1WW7pwUYW5aab8hpX+7SZgK4IsalAPA2skrZb0EuBmYGdTmaeANwFIehWwFngig31Xiptmyi2vM4CYjUrqGlREnJR0K7AbGAPuiIgDkm6pv34b8CfAZyU9ykKT4Ici4sdp911k/dwpOzuv3DyZrtnZMhmoGxG7gF1N225reHwU+JUs9lUGaQZKummmvPrtY/SsEFZWnotvBJyNZ0n6mUzXzYJWZg5QI+BsPEvSTx+jb3aszDwX3whcMDHOsePzS7Y7G6/a+ulj9M2OlZkD1JDt2DfL8yeWTqIxvkzOxrOe+xg99MDKzE18Q7Zt9yHmTy0d3nX+uee4Y9t65qEHVmauQQ1Zq6aXYy8sbfIz66TbZkFn+lkROUANWRZNMr7YWKNOzYJe/8mKyk18Q5a2ScZpxdYrZ/pZUTlADVnaiTp9sbFeOdPPispNfCOQZjYIX2ysV870s6JyDapg+pltwKrNmX5WVA5QBeOLjfXK6z9ZUbmJr2A8o7n1w5MMWxE5QBWQLzZmVgVu4jMzs1xygDIzs1xyE59ZBSTNPgLuy7R8c4AyK7mkqY623LMfAuZPx5ltnv7I8sZNfGYllzT7yPypOBOcFnlGEsubTAKUpBskHZJ0WNLWFmWulfSIpAOS/nsW+zWzznqZZcQzkliepA5QksaATwJvAS4H3inp8qYyk8CngLdFxGuBm9Lu18y608ssI56RxPIkixrUVcDhiHgiIk4AdwEbm8q8C7gvIp4CiIhnM9ivmXUhafaR8TExvkxnbfOMJJY3WQSoKeDphudH6tsa/TPgQknflLRX0ntavZmkzZJmJM3Mzc1lcHhm1ZY01dG2d1zJtpuu9PRHlmtZZPEpYVvzmubnAG8A3gRMAP9L0kMR8YMlvxixHdgOUKvVlq6NXhBeVNDypNXsI/4/aXmWRYA6AqxseH4JcDShzI8j4nngeUnfAq4ElgSoImoORm+8bDn37p31CqZmZilk0cT3MLBG0mpJLwFuBnY2lfkS8C8lnSPpPOBq4GAG+x65pBVu73zoKS8qaGaWUuoaVESclHQrsBsYA+6IiAOSbqm/fltEHJT0NeB7wGng9oj4ftp950HSGJNW7ZJO4bUycPO1DUsmM0lExC5gV9O225qebwO2ZbG/POkl6DiF14ouaVYKN1/boHgmiZRaBZ3mzBGn8FoZJLUYuPnaBsUBKqVWK9y++5pVTuG10mnVYuDmaxsETxabkle4tSpZMTnBbEIwcvO1DYIDVAa8wq1VxZYNa8/qgwI3X9vgOECZWdfcYmDDpIj8TtZQq9ViZmZm1IdhZi045dz6IWlvRNQ6lXMNysz64pRzGzRn8ZlZX5xyboPmAGVmfXHKuQ2aA5SZ9aVVarlTzi0rDlBm1pdWg9Sdcm5ZcZKEFYqzxvLDKec2aA5QVhjOGhuOXm4CPEjdBskBygqjVdbYB+/eDzhIZcE3AdZslK0W7oOywmiVHXYqgg/f9yg79s0O+YjKx6nj1ihpQdZhftccoFrYsW+W9dN7WL31ftZP7/HFLwfaZYf5IpoNp45bo1HfsLiJL0FSM8eWe/bzkZ0HeO74/JlqLriDeJiSJipt5Itoep6t3BqN+obFNagESXcN86eCY8fnz1Rzt/zVfrbcs39kVd8q2rRuio+9/QrG1Lwc5AJfRNNz6rg1avWdmjxvfCj7d4BK0M3dwfzpYP7U2RPtuplpMBqbW7ftPsQ7r17pi+iALN4EeLFNg4UblvGxpTeE//j/Tg7lZjyTACXpBkmHJB2WtLVNuV+UdErSO7LY76CkuRN3M1O2kjpp7907y2+8YcoX0QHZtG6Kb2+9jien38q3t17nv2uFbVo3xUtfsrQnaP50DOVmPHUflKQx4JPA9cAR4GFJOyPisYRyfwbsTrvPQevU19GOm5my1aqT9sHH5/j21utGdFRm1fHc8fnE7cO4Gc+iBnUVcDginoiIE8BdwMaEcr8H3As8m8E+B6q5mePC88YZX3Z2NXd8mZZUfd3MlL1Rd9KaVd0o51zMIotvCni64fkR4OrGApKmgF8HrgN+sd2bSdoMbAZYtWpVBofXn+YR8kmD1cBZfIPmrDKz0UpqURrWzXgWASoppap5md6/AD4UEafUIgPrzC9GbAe2w8KKuhkcX9+6GUHtgDRYo/xymNlo51zMIkAdAVY2PL8EONpUpgbcVQ9OFwE3SjoZETsy2P9AeMqX0Vu8QTg+f4oxiVMRTLmmOhKepLfaRjXnYhYB6mFgjaTVwCxwM/CuxgIRsXrxsaTPAl/Jc3CC9iOo/cUcvOYbhFMRZ2pO/vsPl2/WyqkINx2pkyQi4iRwKwvZeQeBuyPigKRbJN2S9v1HxZ3zozXqKVbsRT4X5TPqOfa6lclURxGxC9jVtO22FmXfl8U+B82d86PlG4T88Lkon6K0EHkmiRY85ctoeTnx/PC5KJ9ONx15mSzbAaoFT/kyWr5ByA+fi/Jpd9ORp+Y/RYw0k7utWq0WMzMzoz4MG5EidOJWhc9FuTQnvsDCTcfH3n4F23YfSuzemJqcyGz2Fkl7I6LWsZwDlJlZ9bS66Vi99f4lA1lhYcDrk9NvzWTf3QYorwdlZlZBrcY25SlBzH1QZmZ2Rp76HCtVg3I7ulm2/J0qn1FObdSsMn1Q7ToF/YUy613Sd0rAu69ZxZ9uumJ0B2a55z6oJkUZmGZWFEnfqQDufOgpAB58fG7kd+BWbJUJUB4NXx5uVsqHVt+dxSC12DbjufusX5VJkvBo+HLI0yDCqmv33WnuOPDcffmSl5kiOqlMgMpTZor1zxOX5seWDWsTF4Nrxa0V+VCkm7zKBChPXVQObqrNj03rpnj3NauWBKlWQcutFflQpJu8yvRBQfLANPdnFEueBhEa/OmmK6j9/CvO+g698bLl3Lt31qsg51SRbvIqFaCaeSG24vES8PmTdOPXHLR845cfF0yMc+z4/JLtebzJq3SAcup58eRpEKG1Nqolwq29Hftmef7EySXbx5cplzd5lQ5QRarq2ot88TPrz7bdh5g/tXRyhvPPPSeX36lKByj3Z5hZlbS6+f7ZC/Osn96Tu1aJTLL4JN0g6ZCkw5K2Jrz+bknfq//7G0lXZrHftJx6bmZV0urmW5DLtPPUAUrSGPBJ4C3A5cA7JV3eVOxJ4F9FxC8AfwJsT7vfLDj13MyqJOmmXOR3YHUWTXxXAYcj4gkASXcBG4HHFgtExN80lH8IuCSD/WbCqedmVhVJSUZJ3RyQj774LALUFPB0w/MjwNVtyv828NUM9jsQTj03szJrvilfP70nt33xWfRBJQ0cT1zDQ9IbWQhQH2r5ZtJmSTOSZubm5lIdWD/zTRVplLWZWVp57ovPogZ1BFjZ8PwS4GhzIUm/ANwOvCUiftLqzSJiO/U+qlqt1vdiVf3WhJx6bmZVkuexhVkEqIeBNZJWA7PAzcC7GgtIWgXcB/xWRPwgg3121O8gXKeem1nV5HVsYeomvog4CdwK7AYOAndHxAFJt0i6pV7sj4CfAz4l6RFJ2SyT20a/NaE8V3fNzKokk4G6EbEL2NW07baGx78D/E4W++pWvzWhPFd37UWtMi2dgWlWHorou5tn4Gq1WszM9FfZau6DgoWakMc5FV+rc/sbb5hKnEXb59wsXyTtjYhap3KlXQ/Kg3DLq1X/4he/87QzMM1KpNRz8eW148/SadWPeKpFa4AzMM2KqbQ1KCuvVv2IY0pey9UZmGbF5ABlhdMq0/KdV690BqZZiZS6ic/KqV2mpVdyNSuP0mbxmdloONXfOuk2i881KDPLjCdbtiy5D8rMMuPJli1LDlBmlokd+2ZzvbaQFY8DlJmltti014pT/a0fDlBmllpS094ip/pbvxygzCy1dk14nmLM+uUsPjNLrdXqAVOTE10FJ6emWxLXoMwstTTrqC32X80eO07wYmr6jn2zAzpaKwoHKDNLLc3qAU5Nt1bcxGdmmeh39YB+V7+28nMNysxGqlUKulPTrbI1KHfKmuXDlg1rE1dIdmp6Z2W/jmUSoCTdAHwcGANuj4jpptdVf/1G4AXgfRHx3Sz23Q/PF2aWH+1mp7fW0lzHihLYUgcoSWPAJ4HrgSPAw5J2RsRjDcXeAqyp/7sa+HT950i065TN40my3jV/Ad942XIefHwu91/IqvLq173r9zpWpBv0LPqgrgIOR8QTEXECuAvY2FRmI/D5WPAQMCnp4gz23ZdWna+zx46zfnqP01sLLilt+QsPPeU0ZiuVfpNLipQ1mUWAmgKebnh+pL6t1zIASNosaUbSzNzcXAaHt1S7zldfvIqv3bQ7i/L6hTTrVr/JJUXKmswiQClhW/MqiN2UWdgYsT0iahFRW758eeqDS5I0qLCRL17F1u0XLY9fSLNu9Ts4ukhZk1kEqCPAyobnlwBH+ygzNI2DClvxxau4uv2i5fELadatfgdHp5n1Y9iyyOJ7GFgjaTUwC9wMvKupzE7gVkl3sZAc8VxEPJPBvvu22Cm7fnpP4hxivngVV1LacrO8fiHNepGUXNIpQ69IWZOpA1REnJR0K7CbhTTzOyLigKRb6q/fBuxiIcX8MAtp5u9Pu9+seAxG+SR9AZ3FZ1XQbYZeUbImFZHYFZQLtVotZmZmBr6foowJMDNrp1WL0NTkBN/eet0IjiiZpL0RUetUrrIzSTQqyt2EmVk7RcrQ64bn4jMzK4kiZeh1wwHKSmHHvlnWT+9h9db7PdjaKqtIGXrdcBOfFV6Rpm4xG6QiZeh1wwHKCs9zK5q9qEx96m7is8IrW8ewmS1wDcoKb8XkRKaDrT3swCwfXIOywsuyYzhpJnRPHmw2Gg5QVnj9zkmWpEhLEZiVnZv4rBSy6hh2f5ZZfrgGZdagbAMdzYrMAcqsQdkGOpoVmZv4zBqUbaCjWZE5QJk1KdNARyuPKg5/cIAya1DFi4DlX1Wn83IflFmdx0BZXlV1+IMDlFldVS8Cln9VHf7gJj6zum4uAm4CtFHIejqvonANyqyu0xioHftm2XLP/rOaALfcs99NgDZwVR3+kCpASXqFpK9L+mH954UJZVZKelDSQUkHJH0gzT7NBqXTReCjXz7A/Kk46/X5U8FHv3xgaMdo1dTvdF5FX8gzbQ1qK/BARKwBHqg/b3YS+GBEvAa4BvhdSZen3K9Z5jpdBH72wnzi77XabjZKZUj6SdsHtRG4tv74c8A3gQ81FoiIZ4Bn6o//QdJBYAp4LOW+zTLnMVCWR/2kmZdhIc+0NahX1QPQYiB6ZbvCki4F1gHfaVNms6QZSTNzc3MpD88sO5MT4z1tN8tKPxmmZcj86xigJH1D0vcT/m3sZUeSzgfuBX4/Iv6+VbmI2B4RtYioLV++vJddmA3UR972WsaX6axt48vER9722hEdkVVFP8GmDBMfd2zii4g3t3pN0o8kXRwRz0i6GHi2RblxFoLTnRFxX99HazZCnqfPhqlxSMMyiVMRS8q0CzZbNqw9q1kQipf5l7YPaifwXmC6/vNLzQUkCfgMcDAi/jzl/sxGyn1UNgzNfU5JwalTsCnDDZUi4YN3/cvSzwF3A6uAp4CbIuKnklYAt0fEjZL+BfA/gEeB0/Vf/YOI2NXp/Wu1WszMzPR9fGbD4gG8lqX103sSB+aOSZyOKPz/MUl7I6LWqVyqGlRE/AR4U8L2o8CN9cf/E1BzGbOyqOpEnjY4rfqWTkfw5PRbh3w0o+OZJMxS8hx+lrUyJDhkwQHKLKUypPNavlR1aqNmDlBmKflu17LW79RGZePZzM1SKkM6r+WPM0YdoMxSK0M6r6W3Y98sH/3ygTNzM05OjPORt73W/w9ScIAyy0Cvd7tOSy+XxaVYGme7P3Z8ni1/tR9wNme/HKDMhqAxIF0wMc7zJ06euZg5Lb34tu0+tGQpFoD501GoyVnzxkkSZgPWvOzBsePzSy5mTksvtnYZm87m7J8DlNmAJY2TSuILWXG1y9h0Nmf/HKDMBqzbwOMLWXFt2bCW8bGlE+aML5OzOVNwgDIbsG4Cj9PSi23Tuim2veNKLjzvxbXBJifG2XbTle5/SsFJEmYDljROanyZOP/cczj2wryz+ErC45ay5wBlNmAeJ2XWHwcosyHw3bVZ79wHZWZmueQalJnZgHjGkHQcoMzMyD6YeCHL9NzEZ2aV1zzbx2Iw2bFvtu/39EKW6aUKUJJeIenrkn5Y/3lhm7JjkvZJ+kqafZqZZW0QwcQLWaaXtga1FXggItYAD9Sft/IB4GDK/ZlVyo59s6yf3sPqrfezfnpPV3f0/fxO1Q0imHghy/TSBqiNwOfqjz8HbEoqJOkS4K3A7Sn3Z1YZ/TQ7DaKpqgoGEUy8bHt6aQPUqyLiGYD6z1e2KPcXwL8DTqfcn1ll9NPs5H6P/gwimHjZ9vQ6ZvFJ+gbw6oSX/rCbHUj6VeDZiNgr6douym8GNgOsWrWqm12YlVI/zU7u9+jPoGb78ADtdDoGqIh4c6vXJP1I0sUR8Yyki4FnE4qtB94m6UbgXODlkr4QEb/ZYn/bge0AtVpt6QpgZiXRKa15xeQEswmBpdPSDr3+ThEMYzyRV0XOn7RNfDuB99Yfvxf4UnOBiPhwRFwSEZcCNwN7WgUns6ropq+on2anMvZ75LFfLY/HVEZpA9Q0cL2kHwLX158jaYWkXWkPzqysuukr6qcPo4z9HnnsV8vjMZVRqpkkIuInwJsSth8FbkzY/k3gm2n2aVYG3fYV9dOHUbZ+jzz2q+XxmMrIM0mYjYDHyHQvj3+rPB5TGTlAmY1AGfuKBiWPf6s8HlMZebJYsxFIk9acVfZYUbLQ8rjgYx6PqYwUkd9M7lqtFjMzM6M+DLPcaJ4hGxbu3HtNhEh6Hy9Db8MiaW9E1DqVcw3KrEDaZY/1EkyS3mf+dPCzF+aB/C4NUZRan2XDfVBmBZJV9lg35fOWNu2xR9XjAGVWIFllj3VbPm3adJYzq3vsUfU4QJkVSFbZY0nvkyRN2nTWNR6PPaoe90GZFUhW2WPN73PBxDjPnzjJ/KkXk6bSpE3v2DfLB+/ez6mmJKx++ssWtZpncPK88cT9u6+q+JzFZ2ZAtunrzRmCjQQ8Of3Wvt53yz37zwqisJB9uO2mK88ca1aZjjY4zuIzs55kNUVSUl9Ro36bDTetm+IjOw9w7Pj8WdvnT8dZtbKsMh1t9NwHZWZ9aZUA0a5PKO1sC881BadFjftstf/ZY8dTJ2rYcLkGZWY9a25Gaxw31aqvaExK3czWzXpXrco0H6drU/nnGpSZ9axdM1qrTMP/+K+vTB0Uusli7JSh6NT04nANysx61i7le5Dz1HXz3o1lWtWknJpeDA5QZtazTk1tg1yTqpv3XiyzfnpPxyZByy838ZlZz4qy3ERRjtOSuQZlZj0rynITRTlOS+aBumZmNlTdDtRN1cQn6RWSvi7ph/WfF7YoNynpHkmPSzoo6ZfS7NfMzMovbR/UVuCBiFgDPFB/nuTjwNci4jLgSuBgyv2amVnJpQ1QG4HP1R9/DtjUXEDSy4FfBj4DEBEnIuJYyv2amVnJpQ1Qr4qIZwDqP1+ZUOafAnPAf5W0T9Ltkl7a6g0lbZY0I2lmbm4u5eGZmVlRdQxQkr4h6fsJ/zZ2uY9zgNcDn46IdcDztG4KJCK2R0QtImrLly/vchdmZlY2HdPMI+LNrV6T9CNJF0fEM5IuBp5NKHYEOBIR36k/v4c2AcrMzAzSN/HtBN5bf/xe4EvNBSLi/wJPS1ocGfcm4LGU+zUzs5JLG6Cmgesl/RC4vv4cSSsk7Woo93vAnZK+B7wO+Pcp92tmZiWX64G6kuaAvxvQ218E/HhA7z1q/mzF5M9WXGX+fIP4bD8fER2TDHIdoAZJ0kw3I5mLyJ+tmPzZiqvMn2+Un82TxZqZWS45QJmZWS5VOUBtH/UBDJA/WzH5sxVXmT/fyD5bZfugzMws36pcgzIzsxxzgDIzs1yqTICSdJOkA5JOS2qZMinpBkmHJB2WVIgpmXpYl+tvJT0q6RFJuV4JstN50IL/XH/9e5JeP4rj7EcXn+1aSc/Vz9Mjkv5oFMfZD0l3SHpW0vdbvF7k89bpsxXyvElaKenB+lp9ByR9IKHMaM5bRFTiH/AaYC3wTaDWoswY8H9YmIH9JcB+4PJRH3sXn+0/AFvrj7cCf9ai3N8CF436eLv4PB3PA3Aj8FVAwDXAd0Z93Bl+tmuBr4z6WPv8fL/MwuTQ32/xeiHPW5efrZDnDbgYeH398cuAH+Tl+1aZGlREHIyIQx2KXQUcjognIuIEcBcLa17lXcd1uQqmm/OwEfh8LHgImKxPWJx3Rf0/1pWI+Bbw0zZFinreuvlshRQRz0TEd+uP/4GFBWWnmoqN5LxVJkB1aQp4uuH5EZaeqDzqZl0ugAD+m6S9kjYP7eh61815KOq56va4f0nSfklflfTa4RzaUBT1vHWr0OdN0qXAOuA7TS+N5Lx1XG6jSCR9A3h1wkt/GBFLZlpPeouEbbnIw2/32Xp4m/URcVTSK4GvS3q8fleYN92ch9yeqw66Oe7vsjBX2T9KuhHYAawZ+JENR1HPWzcKfd4knQ/cC/x+RPx988sJvzLw81aqABVt1q7q0hFgZcPzS4CjKd8zE+0+W5frchERR+s/n5X01yw0N+UxQHVzHnJ7rjroeNyNF4eI2CXpU5IuiogyTEZa1PPWUZHPm6RxFoLTnRFxX0KRkZw3N/Gd7WFgjaTVkl4C3MzCmld513FdLkl3BktoAAABG0lEQVQvlfSyxcfArwCJ2Ug50M152Am8p55ddA3w3GIzZ851/GySXi1J9cdXsfA9/cnQj3QwinreOirqeasf82eAgxHx5y2KjeS8laoG1Y6kXwf+C7AcuF/SIxGxQdIK4PaIuDEiTkq6FdjNQrbVHRFxYISH3a1p4G5Jvw08BdwEC+tyUf9swKuAv65/f84B/jIivjai422r1XmQdEv99duAXSxkFh0GXgDeP6rj7UWXn+0dwL+VdBI4Dtwc9VSqvJP0RRay2S6SdAT4Y2Acin3eoKvPVtTzth74LeBRSY/Ut/0BsApGe9481ZGZmeWSm/jMzCyXHKDMzCyXHKDMzCyXHKDMzCyXHKDMzCyXHKDMzCyXHKDMzCyX/j9o/mmu73//HgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[:, 0], X[:, 1])\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## K Means and Agglomerative Clustering\n",
    "\n",
    "There exist several algorithms for partitioning data into partitions, two of the more common of which are called K Means and Agglomerative Clustering.\n",
    "\n",
    "The K Means algorithm approaches the clustering problem by partitioning a set of data points into disjoint clusters, where each cluster is described by the mean of the samples in the cluster. The mean of the samples in a particular cluster is called a centroid; the K Means algorithm finds the centroids and associates data points with centroids in such a way as to minimize the within-cluster sum-of-squares.\n",
    "\n",
    "For more information on the K Means algorithm and its implementatin in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#k-means\n",
    "\n",
    "In the code cell below, we instantiate the `KMeans` algorithm from the `sklearn.cluster` module and apply it to our data using the `fit_predict` method. We see that `KMeans` identifies two centroids; one located at about (-0.23, 0.56) and the other located at (1.17, -0.05)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.17408322 -0.05027964]\n",
      " [-0.23011109  0.56790752]]\n"
     ]
    }
   ],
   "source": [
    "from sklearn.cluster import KMeans\n",
    "\n",
    "km = KMeans(n_clusters=2, random_state=0)\n",
    "y_km = km.fit_predict(X)\n",
    "print(km.cluster_centers_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Agglomerative Clustering algorithm behaves a little bit differently and does not identify clusters using centroids. Instead, it recursively merges the pair of clusters that minimally increases a given linkage distance. Put another way, the Agglomerative Clustering algorithm identifies the two data points that are \"closest\" out of all the data samples. It then takes those two data points and identifies a third data point that is \"closest\" to those two data points. The algorithm continues in this fashion for each data point; finding the next data point that is \"closest\" to the preceeding cluster of data points, where the definition of \"closest\" depends on the distance metric chosen.\n",
    "\n",
    "For more information on the Agglomerative Clustering algorithm and its implementatin in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#hierarchical-clustering\n",
    "\n",
    "Below, we instantiate the `AgglomerativeClustering` algorithm from the `sklearn.cluster` module and apply it to our data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cluster import AgglomerativeClustering\n",
    "\n",
    "ac = AgglomerativeClustering(n_clusters=2,\n",
    "                             affinity='euclidean',\n",
    "                             linkage='complete')\n",
    "y_ac = ac.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can visualize the results of both algorithms applied to the data. Visually, we see that neither algorithm ideally clusters our data. The ideal algorithm for this unique set of data would recognize that both sets of samples are generated from two different equations describing two different half circles."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAADQCAYAAAAK/RswAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsnXl8FPX5+N8PBARzbjgNV0BFQQuIgKJUDg/AA/uttD9PtEApRc4EIR5VpGpBEQU8kAZaj289vmJbb9AiINYLC6ggl4gQI3dCEoSawPP7Y2bjZjO72U12s9fn/XrNK7szn8/MZyazzzzzfJ5DVBWDwWAwGAyGeKJBpAdgMBgMBoPBEGqMgmMwGAwGgyHuMAqOwWAwGAyGuMMoOAaDwWAwGOIOo+AYDAaDwWCIO4yCYzAYDAaDIe4wCo4hKhGRnSJySRSM4wYRWR7pcRgM9YmI3CIiayI9jtoQqd+siAwQkYL6Pq4TIvKWiNwc6XFEGqPgRDHeD3kRuVZEikSkv0PbASKiIvKK1/ru9vqV9TDkgBGRNBF5VER2iUiZiGy3vzcP4TFmiMhzddmHqv6vql4WqjEZDKFGRFbacuGkSI+lvhGRbFu+JbnXhfM3KyJ9RORNESkWkUMi8omI/CbEx6jzy52qDlXVp0M1pljFKDgxgq2NPw5coaqrfDTbD1wgIs081t0MbA33+IJBRBoD/wLOAoYAacAFwEGgTwSHVgVPoWkwRCMikg38HFBgWEQHEwZEpGGkx+BGRPoCK4BVwGlAM+D3wNBIjssTsTDPdRtzIWIAERkDPAwMVtV/+2n6I/AP4Fq7X0Pg18D/eu3vTBF5x34D2SIiv/bYdoWIrBOREhHZLSIzPLa535Zuti0vB0TkTo/tfURkrd13r4jM9THOEUB74H9UdZOqnlDVfar6R1V90+H8/yoi93l8r2IKFpHpIvKdiJTa53OxiAwB7gD+n20h2mC3TReRxSLyvd3nPrcQtc3yH4jIIyJyCJjhbaq3z3+siGyz35ofFxFxX28Redi+Lt+IyHjvt0uDIcSMAD4C/or1MlOJiDQTkdfs3+On9r3ueS9fZv9eDovIEyKySkRGOx1ERC6w93HY/nuBx7aV9r7/bf/WXrOP/b8ex872aO9P/vxVRJ60rSRHgIH+ZBKw2v5bbB+7r+dvVkQWisgcr3P5p4jk2J+zRGSpiOy3f7MT/Vzrh4CnVXW2qh5Qi89U9ddOje3f/mle53af/bm5iLwuP1mC3heRBiLyLJZsfM0+n2l2+/Pt61ssIhtEZIDX9b9fRD4AfgA62etG29tvEZE1IjLHllnfiMhQj/4dRWS1LT/ftWVanSzfUYOqmiVKF2AnsBTYC3Svoe0AoADLEvKxve5yYBkwGlhpr0sGdgO/AZKAnsAB4CyP/fwMS/ntZh/7F/a2bKw3xT8DTYHuwH+BLvb2D4Gb7M8pwPk+xvoClqCo6dwvsT//FbjP+1ztz2fY55PlMcZT7c8zgOe89vsP4Cn7OrQEPgF+Z2+7BagAJtjXpqm9bo1HfwVeBzKwBNF+YIi9bSywCWgLuIB37fZJkb6XzBKfC7AdGAecC5QDrTy2vWAvJwNd7d/JGntbc6AE+KV9r0+y+4+2t9/i0TYTKAJustteZ39vZm9faY/jVCDd/g1sBS6x2z8D/MVuW5P8+StwGLgQSwY1ITCZlORx3p5jv8g+ntjfXcBRIMve32fA3UBjoBOwA+tF0vs6nwwcBwb6+V8MwJZL9ncFTvP4/ldsOQb8CVgINLKXn3uMcSe27LO/t8Gybl9uj/lS+3sLj+u/C8sinmTvb6XX/7Ic+C3QEMvqVOhxvA+BOfY16Id1Xzzn6zxjaTEWnOjnUqw3tC8CaayWhSdTRM7Aert7xqvJlcBOVf2Lqlao6n+wlKjhdv+VqvqFWlaVz4HnAW+fn3tV9aiqbgA2YCk6YP2IThOR5qpapqof+RhmM+D7QM4nAI4DJwFdRaSRqu5U1a+dGopIKyxz8mRVPaKq+4BHsC1eNoWqusC+Nkd9HHOWqhar6i7gPaCHvf7XwDxVLVDVImBWCM7PYHBERPoBHYCXVPUz4GvgentbQ+Aa4B5V/UFVNwGePhmXAxtV9RVVrQDmA3t8HOoKYJuqPmv/Lp4HNgNXebT5i6p+raqHgbeAr1X1XXvf/wecY7fzK39s/qmqH9gy6FiAMskX72MpGj+3vw8HPlTVQqA3lpIwU1V/VNUdWC9v1zrsx4WlXIRKbpUDpwAdVLVcVd9XW9tw4EbgTVV9074G7wBrsf6Hbv6qqhvta1rusI9vVfXPqnoc6z44BWglIu2xrsPd9jVYA7waonOMOEbBiX7GAp2BfPdUSAA8C4wHBgJ/99rWATjPNnUWi0gxcAPQGkBEzhOR92yT7WH7+N6Ov56C8Acsaw3AKHusm22z9JU+xncQ6wdWZ1R1OzAZy1qzT0ReEJEsH807YL3dfO9x7k9hWXLc7A7gsL7OP8urfyD7Mhhqy83AclU9YH//Gz9NU7XAepv3dT9WuVfth6uvCKAs4Fuvdd9iWRbc7PX4fNThu/s34lf+OIwzUJnkiH1eL2BZncBSAN1T9h2ALK+x3AG0cthVEXCCEMktrOmu7cByEdkhInl+2nYAfuU1zn5eY6lJ1lTKLFX9wf6YgvW/PeSxLpB9xQxGwYl+9gEXY72BPBFgn2exzNZvet24YN28q1Q1w2NJUdXf29v/hqXBt1PVdCwzakCKlapuU9XrsBSG2cDLIpLs0PRdYLCPbU4cwTIRu/EUhqjq31TV/Tar9rGxP3uyG2tKrbnHuaep6lmeuwtwTE58jzU95aZdHfZlMPhERJpiWQz7i8geEdkDTAG6i0h3rKnTCnzfj1XuVfvlybOtJ4VYvy1P2gPf1WLoNckfqP4b9CeTAvm9Pg8MF5EOwHlYFiP3WL7xGkuqql7uvQNbjn6IZRULlB/wIbdUtVRVc1W1E5YlLEdELvZxTruBZ73Gmayqnhbi2sqt77Es/p7jjBu5ZRScGMA2pw4ChojIIwG0/wbLhHunw+bXgc4icpOINLKX3iLSxd6eiqXRHxORPtgm70AQkRtFpIWqngCK7dXHHZo+i/WjXWo7HDawnRLvEJFqwgVYD1wuIpki0hrLYuM+5hkiMkisENljWG+L7mPuBbLFjipQ1e+B5cDDYoWpNxCRU8Uh7L6WvARMEpE2IpIBTA/Rfg0Gb36BdZ93xZoi7QF0wZqSGWFPRbyC5Sh/soiciTVl7eYN4Gci8guxnOBvxevFwYM3sWTG9SKSJCL/zz7u67UYd03yxwl/Mmk/lmWlk6/OqrrObpcPLFNVt2z6BCgRK0ihqVhBAmeLSG8fu5oG3CIit4kdqSpWGo4XfLRfD1xv73cIHtNqInKliJxmK5YlWP9LT7nleT7PAVeJyGB7X03ECrTwpZAGjKp+izXdNUNEGosVKXZVDd1iBqPgxAiquhtLyRkuIn8KoP0aWzHyXl8KXIY1z1yIZbqcjeXHApblZ6aIlGI5370UxDCHABtFpAyYB1yrqsccxvBfLAfEzcA7WD/wT7DMzh877PdZLF+fnVgKyose207C8nU5YJ9LSywzM1hz/wAHReQ/9ucRWM50m7DMzi8TOrPzn+3xfQ6sw3owVOCs5BkMdeFmLL+XXaq6x70AjwE32ErLeCyn3z1Yv6HnsSyY2NNavwIexJoy7or1oPuv94FU9SCW70yu3XYacKXH1FjABCB/nPApk2zLyv3AB/b0zfk+9vE8lsz5m0ff41gP8x7AN1gyJB/rmjmN/d9YMngQsEOsSMtFWL9zJybZ+3dPw/3DY9vpWJbsMizL0BOqutLe9ifgLvt8ptqy/2osubYf6+XwNkL3/L4B6Iv1v70PS75Wuw9iEbcXtcFgCDF2KOZCVfU27xsM9Y6IzAZaq2q1DLe2lbMAuEFV36v3wRmiBhF5EdisqvdEeix1xVhwDIYQYZu5L7fN+G2Ae6ju5G0w1Av29G83seiDFQTwd4/tg0Ukw57evQPLr8VX5KMhTrGnCE+1p+yHYFmL/lFTv1jAJCAzGEKHAPdimXiPYvk53B3RERkSmVSsqZksrGCFh4F/emzvizVl456y/YWf1AiG+KU1lr9WMywr3u9tv6WYx0xRGQwGg8FgiDvMFJXBYDAYDIa4IyanqJo3b67Z2dmRHobBYPDis88+O6CqLSI9jlBhZI3BEH0EKmdiUsHJzs5m7dq1kR6GwWDwQkS8M97GNEbWGAzRR6ByxkxRGQwGg8FgiDuMgmMwGAwGgyHuMAqOwWAwGAyGuCMmfXAMBoPBEPuUlZVRWFhIVlYWKSkpNXeIQcrLyykoKODYsWpVaww10KRJE9q2bUujRo1q1d8oOAaDwWAIGYEoLRUVFUybnsfixYtJc2VSUnSIUaNG8eDsWSQlxddjqaCggNTUVLKzs7FqaxoCQVU5ePAgBQUFdOzYsVb7MFNUhmpkpqUhItWWzLS0SA/NYDBEKRUVFeTkTqVN23ZcPHgobdq2Iyd3KhUVFdXaTpuex6pP1jL3tfeY9/Ya5r72Hqs+WcuUnFy2bt1KWVlZBM4gPBw7doxmzZoZ5SZIRIRmzZrVyfJlFJwEIFiFpai0FIVqS1Fpaf0N2mAwxBS+lJZp0/OqtCsrK2Px4sWMe+BRXC1bAeBq2YpxDzzKokVPMfDSwX6Vo1jEKDe1o67XzSg4CUAoFRZj0TEYDN74U1qWLFlSxSJTWFhImiuzsp0bV8tWuFq25vannvWpHBkMwWAUHENQGIuOwWDwxlJaXI5KS1qGi8LCwsp1WVlZlBQdomjf3ipti/bt5UjJYVwtW/tUjmpLWVlZ3E191ZUZM2YwZ86coPsVFxfzxBNP1Pn4jz32GKeddhoiwoEDB+q8PyeMgmMwGGICEVkiIvtE5Esf20VE5ovIdhH5XER6emwbIiJb7G3GLFAHvJWFiooK5i94jL3fFzoqLSXFRWRlZVWuS0lJYdSoUTxxx+TK9kX79rIgbxKDrrmWpsnJgLNyFKyiEoxfUDQRzQpZbRQcVeXEiRNV1l144YW8++67dOjQIZTDq0JIFBwjeKKXaJ1GMo7MhlrwV2CIn+1DgdPtZQzwJICINAQet7d3Ba4Tka5hHWkc4ktZmHrbND5c/zmD/uf/sSBvUhWl5Yk7JjNy5Mhq0VQPzp5F/z69yB02iElDLmTcpefTun02N+bcUdnGUzmqraISqF9QtBAuheyZZ56hW7dudO/enZtuuqna9gEDBlSWJDlw4ADu+msbN26kT58+9OjRg27durFt2zby8vL4+uuv6dGjB7fddhsADz30EL1796Zbt27cc889AOzcuZMuXbowbtw4evbsye7du6sc85xzzqk8TthQ1TovwEVAT+BLH9svB94CBDgf+Nhe3xD4GugENAY2AF1rOt65556rhqq4UlOd3Gw0yf6rDov17w98Xy6H/j7bpqb6HW+wYzLEBsBaDYFM8bUA2X7kzFPAdR7ftwCnAH2BZR7rbwduD+R4iSBrSktLdcuWLVpaWuq33ZScXO3Zr7/mr16nSzcXav7qdXrOhRfpySmpmr96nb705S4dNnKspqRnaKt2HbTxSU10/ISJWl5e7vN47s+3jp9Qbd89+/XXKTm5Po/tud3XeaWlZ1T2cS/5q9dpeoarxvMNFZs2bQq4bW3Osya+/PJL7dy5s+7fv19VVQ8ePKiqqvfcc48+9NBDqqrav39//fTTT1VVdf/+/dqhQwdVVR0/frw+99xzqqr63//+V3/44Qf95ptv9Kyzzqrc/7Jly/S3v/2tnjhxQo8fP65XXHGFrlq1Sr/55hsVEf3www/9jq9Dhw6VY3PC6foFKmdCYsFR1dXAIT9Nrgaescf2EZAhIqcAfYDtqrpDVX8EXrDbGoLElyNxBeDC0iy9F1dqquO+DpWUeD84UJz/wSbiyhBFtAE8XxML7HW+1jsiImNEZK2IrN2/f39YBhppysrK2LRpE+MnTAzIWuDLifjWP82joqKcJskpNExK4uZpd7Nwxafcteg5WrRqzYTxt1bmtXGyTtx9zww6derEo4/MrbToTB7Sj9xhg+jfpxcPzp4VlAOzJ/6cmb2nvqKB2p5nTaxYsYLhw4fTvHlzADIzMwPu27dvXx544AFmz57Nt99+S9OmTau1Wb58OcuXL+ecc86hZ8+ebN68mW3btgHQoUMHzj///FqNOxTUlw9OnQVPIggdN6GevjlEdQUELEUmkGO6UlODUpAMhgjhFFOqftY7oqqLVLWXqvZq0aJFyAYXDXgqGQMuuYz8/Hz6Dh3G3NdX+p2+8acsnJySxs7NGyvXNU1OpmlyCmUlh6v43vibLkpKSmLuw3Mo2L2Ld5e9RcHuXcx9eA7Hjh3jgw8+ICU9I2hFxZ8zs7dfUDQQLoVMVWsMt05KSqr0kfHMO3P99dfz6quv0rRpUwYPHsyKFSsc93/77bezfv161q9fz/bt2xk1ahQAybY/VaSoLwWnzoInnoUOVFUwPK0fLurHKuLPEuNt0XEvngqSwRAFFADtPL63BQr9rE84pk3PY+XHnzL3tfdY+N6nPL783+zZtZPn5j7g11rgT1n479Ej/N9jc/z63gRqnUhJSaFz5840adKkUhEbPfb37N/7PYtm5HHcw7pUtG8vh4sPceTIEUfrhi9nZl9+QZEmXArZxRdfzEsvvcTBgwcBOHSoui0+Ozubzz77DICXX365cv2OHTvo1KkTEydOZNiwYXz++eekpqZS6vEsGjx4cJX/4Xfffce+fftqNdZQU18KjhE8NeBTwYjoqJzxVMbgJ4tO4IZPgyEsvAqMsIMazgcOq+r3wKfA6SLSUUQaA9fabROK4uJiFi5cyK1/mldFyZgwax4rlr7A0SNHfFoL/CkLY8aMYVDf8xynl9wEa52oau35gCff+Yhd27ew+L67Ko/9SO7vqag4zrBfDvc5vebpzOxrbNFCuBSys846izvvvJP+/fvTvXt3cnJyqrWZOnUqTz75JBdccEGVkO0XX3yRs88+mx49erB582ZGjBhBs2bNuPDCCzn77LO57bbbuOyyy7j++uvp27cvP/vZzxg+fHgVBcgX8+fPp23bthQUFNCtWzdGjx5dq/PzSyCOOoEs+Hf+u4KqTsaf2OuTgB1AR35yMj6rpmPFo+Mf/pxuA3DA9dffafF2AvZ7/FqO1d8xa+ucbIhuCKOTMfA88D1QjvVyNAoYC4y1twtWtNTXwBdAL4++lwNb7W13BnrMeJI1t/zmN5rZqnUVh1v30rpDR13w1vt+HXDLy8t1Sk6upme4tF12J03PcOmUnNxKJ2J/zsrBOPz6a9v4pJO0eatTtNFJJ2lKeoY2TUnVYSPH6lMrPvXrjFuTI3Wgjta1IRgn45qucSJSFyfjkFQ1E5HngQFAcxEpAO4BGmE9HRcCb9oCZjvwA/Abe1uFiIwHlmFFVC1R1Y3VDmCoEVdqKuKgNbtSU+t1Ksk955iE5eAM1jTbIagyPjO9ZQgWVb2uhu0K3Opj25tYcighKSsrY+nSVziuStG+vVUsKUX79lJadAiRBn6tBW4/mZn3znAspumeXnLC0zrhnqbyZZ3wZ+1Jz2xGqiuTGU+/zCnZHSvz57zxnDX9lTtsEDPvnVFt/L7GFm1FP2u6xobgCFUU1XWqeoqqNlLVtqq6WFUX2sqN2wRwq6qeqqo/U9W1Hn3fVNXO9rb7QzGeRCRa/GTcpphyj8+1nWYzuXIMhtDgVhradjqdOZPHVJkCmTN5DEkNG3LntVcGNH3jVhaCffAGOl3kzxfl8KGD5Mx9ilOyrerSnlNsTZJTgnbGjdY8ObW9xoaqxFdd+hilpge22yoSzqglfxagSOH2S/LGaZwGg8E3WVlZFB08QOYpbcjK7sTkqwaSnJZO8YF9NGjQkFXvraBLly5hfaAGap3wZe15LG8SJ6emVSo3blwtW5GS4WLn5o1BOeO6HZ/nvvZeNcdnX5YgQ2xhFJwooKi0tDJXjTdJQLn6jGgNGWbKyGCIb/TECSY9uABXy1bcdNsfKNq3B5EG3HbNZbVSbsrKymo1jeJvKsvNg7NnMW16HrnDBpGW4aKkuIibbrqJr7/c4HOK7f8emxOUM24gjs81jdMQ3ZhaVFGCr1w10VgxxVdeHKMtGwzRSWFhIc1atqp8mDdNTiar46mckt2RZi1aBTWtUx/1nZzy4iyYP4/Ro0dXizJ6eMrvQJVBfc8LKjoq1vLkGILHKDhRTrBTRPXht+LL38efeDOJAQ2GyFBWVsaRI0c4XHTQ8WFeGuTDvD79Vrx9Ubz9eHKGDeT87mdTsOtb5j48JyjH4FjLk2MIHqPgRDnBTh1FsnSCv4zHkXB4NhgSGU9Ly7BfDuf48RM8OnVcnR7m4SonECjelp3vdu9myeLFZGRk1Gp/sZQnJ9TMmDGDOXPmBN2vNtXEnbjhhhs444wzOPvssxk5ciTl5eV13qc3RsGJAN5Wlngh1JFcpkSEwVB7vC0tj7y6gkN79zBhaL9aP8yjpb5TqKKMfJWIcLIElZWVsXXr1rArcdFObRQcVa0sBeHmhhtuYPPmzXzxxRccPXqU/Pz8UA4TMApORPC2sgRbDDNRCEZhMsLHYPgJJ0tL86w2zHz2FRolNeKfr7zs92Hui3j1W/GnMNWHz5GbcLkYPPPMM3Tr1o3u3btz0003Vds+YMAA1q61srccOHCA7OxsADZu3EifPn3o0aMH3bp1Y9u2beTl5fH111/To0cPbrvtNgAeeughevfuTbdu3bjnnnsA2LlzJ126dGHcuHH07NmT3bt3Vznm5ZdfXnl+ffr0oaCgoE7n6IRRcKIAt4MxENEcNrFIfQofgyFa8Vbw/SbLc2WSnJxcK+tHIvqt1KfPUThcDDZu3Mj999/PihUr2LBhA/PmzQu478KFC5k0aRLr169n7dq1tG3bllmzZnHqqaeyfv16HnroIZYvX862bdv45JNPWL9+PZ999hmrV68GYMuWLYwYMYJ169bRoUMHx2OUl5fz7LPPMmTIkFqfoy+MgmOIaaI1UZfBUB/4UvBbtmwZNktLIvmtRNrnKBSsWLGC4cOH07x5cwAyMwOvGti3b18eeOABZs+ezbfffkvTpk2rtVm+fDnLly/nnHPOoWfPnmzevJlt27YB0KFDB84//3y/xxg3bhwXXXQRP//5z4M4q8AwCk4Y8WVuDGc4dSL5rcSD8DEY6oIvBX/mH+8Lm6UlGL+VWCdafI7qgqrW6OuZlJRU6SNz7NixyvXXX389r776Kk2bNmXw4MGsWLHCcf+3334769evZ/369Wzfvp1Ro0YBkJyc7Pe49957L/v372fu3LnBnlZAGAUnjPgyN4Zz8iRaSjbUB/EgfAyG2lKTgn/3H+4Kq6Ul1soJ1MZPLx58ji6++GJeeuklDh48CMChQ4eqtcnOzuazzz4D4OWXX65cv2PHDjp16sTEiRMZNmwYn3/+OampqVWqhQ8ePLjKC+V3333Hvn37ahxXfn4+y5Yt4/nnn6dBg/CoIkbBiRCJYmUJJ/EgfAyGQHB6ONek4O/bty9hLC3+qIufXjz4HJ111lnceeed9O/fn+7du5OTk1OtzdSpU3nyySe54IILOHDgQOX6F198kbPPPpsePXqwefNmRowYQbNmzbjwwgs5++yzue2227jsssu4/vrr6du3Lz/72c8YPnx4FQXIF2PHjmXv3r307duXHj16MHPmzJCeN4Dj2360L+eee2618unRCKDqsFBZf9TgSk11MnKpKzU1oP5TcnK1Z7/+mr96nS7dXKj5q9dpz379dUpObphHbnACWKth/O0DQ4AtwHYgz2H7bcB6e/kSOA5k2tt2Al/Y2wIaZ6RlTXl5uU7JydW09Axtm91J09IzdEpOrpaXl2tpaammpWdU3vvuJX/1Ok3PcGlpaWlExx4tBCojSktLdcuWLdWum/t/kJ7h0nbZnTQ9w1X5PwiETZs2BTzWusrDeMTp+gX6+zVCp5YEciMaBadmarpGNV3nugofQ2gJp4IDNAS+BjoBjYENQFc/7a8CVnh83wk0D+aYkZY1NT2cjYLvn0CUQH9KpPe+nBSgmghGwTFUpy4Kjlhta4+INAS2ApcCBcCnwHWquslH+6uAKao6yP6+E+ilqgec2jvRq1cvdcfsRwoRca50DW5hGlCbRKemaxToNaxt4T9DaBGRz1S1V5j23ReYoaqD7e+3A6jqn3y0/xvwnqr+2f6+kxiSNWVlZbRp244HXnwD1RO4WramaXIyRfv2kjtsEAW7d9GkSROmTc9jyZIllUUpR44cyYOzZyXcVJQTW7du5eLBQ5n39ppq2yYP6ce7y95i4VOLWPXJ2iqVy5+4YzL9+/Ri7sPBZ/r15quvvqJLly513k+i4nT9ApUzofDB6QNsV9Udqvoj8AJwtZ/21wHPh+C4UU8iRTRFmlhzeDTUijaAZ7awAntdNUTkZCzL8lKP1QosF5HPRGSMr4OIyBgRWSsia/fv3x+CYdeOXbt20aBRI6b/+nLuHzuC3w3sxdMPziQts1mlE30iRTTVhpr89NLS0mqMxAxFElHzQls76nrdQqHgJJTQcZOEs6Owp1hJpIimSGCyFyccTrGuviTgVcAHquoZMnKhqvYEhgK3ishFTh1VdZGq9lLVXi1atKjbiOvAkwufonX7bOa9vpLHl33AvNdX8u3mTSy+765qTvQpKSlkZWVRWFhofg8e1OQkXFJS4tNROzXDxYSJE+ucRLRJkyYcPHjQKDlBoqocPHiQJk2a1HofoVDzQyF0CkWkJfCOiGxW1dXVdqi6CFgEltm4roOuKxU4n2T8VJaKftq0bUeaK5OSokOMGjXKmOXjnwKgncf3toCvXADX4mUpVtVC++8+Efk7lvW5mqyJBsrKynjmmWeY+9p7VSwLE2bN49bBFzB61OhKa2VFRQXTpuexePFi83tw4MHZs5g2PY/cYYOqTeMdO3as0sLjqeQU7dvLwf17Wbd5W+X/wK0YTZueF9TUVdu2bSkoKCAaXsxjjSZNmtC2bdta9w+FD07A8+K2UPk/Vf2bj33NAMpU1e/dEys+OJlpaY4ptl2pqcal0PAHAAAgAElEQVSKY1PTNfJ3nfNXrwvLnLmh9oTZBycJy9/vYuA7LH+/61V1o1e7dOAboJ2qHrHXJQMNVLXU/vwOMFNV3/Z3zEjJGn++I2MH9mblu8vp2rUrADm5U8PqQxIv+PLTc7p+j98+iS0b1jH/zdXVFB+3/5OZDo8c9emD8ylwuoh0FJHGWG9OrzoMKB3oD/zTY12yiKS6PwOXYUVZxQXhqCsSb9Q0jefLjynl5JNN9uIEQ1UrgPHAMuAr4CVV3SgiY0VkrEfT/wGWu5Ubm1bAGhHZAHwCvFGTchNJ/PmO/Hj0B9q3bw+YbN7B4MtPz6n0RI8zTsPVrLlJIhrj1FnBSSShY6h/vBWgLVu20Da7E0//Z3uVdkbwJAaq+qaqdlbVU1X1fnvdQlVd6NHmr6p6rVe/Hara3V7OcveNVgJNMGeyedcdJ0ft+fPmUVpcZJKIxjghmaBV1TeBN73WLfT6/lfgr17rdgDdQzEGQ2Lg+WbrbTo2gscQT/jzHXFjfg+hw23hceNWML2n/mIlg7HBlGqoNSYEPDLEQ+p0gyEQAgkBN7+H8PHg7Fn07dGNyVf0Z+JlF8R11fR4xbjY1xLjJBw53G+2OVcN5OSUFH4oK6uMGjEY4g1vy4I3gVh6DMHhjkx79tlnyWjWjOJDB7nl5lu4+w93sWPHDpNQNEYwCk4YcaWmIj4ihAx1R1EaJiWhdpxVWVkZ+/btM8LHkFC4LT0z751hsnmHiGnT81j1ydoqIeKPTh1Hm3btyWzewoTixwhmiiqMmER/4cEtfB55bSXzl/2bR15byTtr/k2bdu3rlJDLYIhlTDbv0OArMm3ynCeQBg2YtXQZc197j1WfrGXa9LwIj9bgD6PgGGIKI3wMBkM48ReZlurKpGjfHhOKHyMYBccQUxjhY4h3TAmSyOIvB1FZcRGulq0BE4ofCxgFxxBTGOFjiEfKysrYtGkT4yfUvfaRoW74ikybP30ig665lqbJyZXrTCh+dGMUnADJTEtDRKotmWlpkR5aQmGEjyGeqKioICd3Km3atmPAJZeRn59P36HDmPv6SjPVGkG8sxtPGNqPA4UFXHHjKMB3KL6xvkUZTk6w0b6ce+65Wh+4UlOdKi2oC1TtxbqEhvqkvLxcp+TkanqGS9tld9Kmycma1aGjPrXiU126uVDzV6/Tnv3665Sc3Cr9SktLdcuWLVpaWhqhkcc/wFqNAhkRqiXcsmZKTq727Ndf81evq7x3u19wkQ4bObbye3qGy9yzEcItM4qKiqrInPQMl07JydXy8nJV/UkmpaVnaNvsTpqWnlFluyG0BCpn6lxsMxLUVwE8vwU1PT/H4DWMB9zF81q2bMnMP97HkiVLquUBSUpKMtWW65FwFtuMBOGUNWVlZbRp265KxXCwrAOTrxrIwhWf0jQ5mclD+vHusrf85sIx1A/BFOw0BU/DR30W2zQYIoI7LDYjI8NvxlfPnBbz3l5jTP+GqMCfw3xKhouifXvMVGuU4RSKbwqeRi9GwTHEDUb4GGKJmhzmRRqYkgsxgCl4Gr0YBccQ1xjhY4hWfDnMz5k8hkaNGnHntVea2kcxgD9F9XDxIWN9iyAhUXBEZIiIbBGR7SJSze4vIgNE5LCIrLeXuwPtG62Y4pqxgT/hY0z/sUe8yRrvaJ3cYYO45MK+rFm10rG4piH68KmoTvkdFRXHufueGSbUP1IE4onsbwEaAl8DnYDGwAagq1ebAcDrtenrtEQ8iio1tV6ObwgNTpEqZ/bsrb36nGeiHEIMYYyiimdZYyL8Yht3FFXT5GTNbNVak9PSddjIsfrUik8dIzoNdSNQORMKC04fYLuq7lDVH4EXgKvroW/YMbWk4oMHZ8+iwY/HuHXwBYy7tC+TrxpIh85dON7oJONoHFvErawxdaRim6SkJGbeO4OkhknkzF3IU++t5eZpd9M8q43x94sgobB9tgF2e3wvAM5zaNdXRDYAhcBUVd0YRF+DodYcO3aMrVu28PDf30X1BK6WrWmanEzRvr3kDhvEzHtnmAdLbFAvskZExgBjANq3bx+CYSc2vkKr4wHPcyssLCQ9sxldzu1TpY2nv58J9a9fQmHBEYd13olh/gN0UNXuwALgH0H0tRqKjBGRtSKydv/+/bUerMEi3jMze2YUdTsan5LdkayOp1ZmOzaOxjFHvcgaVV2kqr1UtVeLFi1qPdhExzNL88BLB9P6lCzGT5gYF/4onufmLqsxf8FjlBQdNP5+UUQoFJwCoJ3H97ZYb06VqGqJqpbZn98EGolI80D6euzDCJ0QUlRaWt25yF4fyxjBExqiNOV8vciacBKl1zUsTJuex8qPP6Xv0GEcLi4i1ZVplaK4sF/MKzlOubU+XP85nc84s5qzsQn1r059vWCHQsH5FDhdRDqKSGPgWuBVzwYi0lpExP7cxz7uwUD6GsJHJj9Fg7kXIKatOEbw1A0nBTGKCj7GrKzxvq5ZbdsyctQoiouL41Lpceefatkumz27djLv9ZU8/s6HPL7835RVnGBKTm6kh1hr/OXW2rplCxec071KVJwJ9a9Ofb1g19kHR1UrRGQ8sAwrUmGJqm4UkbH29oXAcOD3IlIBHAWutT2hHfvWdUzBkJmW5nhRXampce9MXISzjV5i1IrjFjyeqe/dgifnqoGMGDGC3GGDqpVzMPyEp4LoTjn/+O2TmDY9L+Ip52NZ1jhd14cmjyGrbTtQpVnLVnFVQqSwsJCTU9NY/dorPLR0WZXf49RHFzHlqgH86YH7Y/Llwl9urXRXJhPG38qfHrg/bv2OYomEr0XlVG8qE+vh7008KT32S67vWlsxeF9s3bqViwcPZd7ba6ptc9fzcTsDGsFTHX+1kSYM7UdhQQEZGRl+92FqUVXH13VdNCOPb7d+xdRHF8VV/aKKigqm5OTy1KKnSEnL4Mf/HuPi4ddxY84dNLQVt0lDLuRfy96OSadbf7+T3GGDKNi9y8iWGvBb5zGAZ4+pRVUH3JaNePNPgZ/mPt24p6YyIzai0BFIUj8Tjusbf2+mTVNSmZKTE6GRxTZO1/XokSOsefOflcoNxE8JkWnT8/j3ug08+c5H5L+/jnmvr+TbzZt4bu4DgPV7LC0ujlnfN1+J/cyUd/RhFJwEw+fcZ0RHFRqM4KkbWVlZHPbhjP3fo0d55ZW/x/SDN1I4Kd5F+/aQmp4RdyVEfPmnTJg1jxVLX+D7nd/Exe/RKQO18bWJPoyCY4grjOCpPSkpKQy/5hrmTPldFQVxQd4kLh5+HemuzJh98EYSt+L9+O2TKq+rSAMO7d8Xd5F9/qyAjRqfRN6vhsbF7zEpKYm5D8+hYPcu3l32FgW7dzHz3hns2LHDvAQEQBLVA1yE0CTm88T44DjMBQrx5Zviib+5Tzfx4GvkmYALMH43AVJcXEybdu2RBg1IdWVSVlzEoGuu5YobRzHtl5fV6F9gfHCcOXbsGD/vP4ANG9aTnJbOkZLDpDdrQVaHjkycPT9ufHCKi4vJatOWBW+vqeafMuWqAWzdvJnWrVtHcIShp6KigmnT81i8eDFprsy4chYPF/Xlg5PwV9+tSSYCNYV/x7ry5klKSgqdOnUygicAPJXBjIwMfjdmDCs+/JhfjZ9K9plncexIWVxMK0SSO+68ixONm/DIP1egegKRBjx2x2S2fb6OCUP70axFK0rjILJv5h/vw9WiJfOnT6yiuD06dRyjR42OO+UGnCPknrhjclREHkYaX1HK9SV9jQUngaKo4jFyyh85uVNZ9cnaSn+AeHhDDiW+3jwfuP8+7rjzLpYsWVItpL4mxdBYcKrjL+omZ9hAtnz1FSUlJTFvYXSf50N/f4c3nlvMiqUvkJLhorToEKhSsOvbGqPwYg0TUeWfQGYMPAn0GWssOHXgEPH5wAdw4XxzuVJT63soYcVfThxTf8rC15vnbdOmM2H8rUy7bWpcPHgjjd+8KRmZlJSUxGS4tDfu82ye1Yabp93Nr2/NpWjfHlwtW3P7NYPZt29f3Ck4/v63pv5U4ITLeGCcjBOMQ1SPoAJi2jLlRCCCJ5Hxl4110aKnGHjpYM44swsLn1pEkyZNIjza2CaQ9AXxgPd5Nk1OJqvjqRw7UhZX5+lJovxvw0F9pGAxCo4hLjGCxz/+FEBXy9bc/tSzzH3tPVZ9spZp0/MiNMrYxLv0QqKkL0iU8/QkEc85lkh4BceVmuoYrhZvUzaJhhE8/vGnAB4pOYyrZeu4STxXX/ir45Uo6QsS5Tw9cTrnvj26cf1117Ju3Trz24kgCeFknMj1pjxJtOvgdqL1dJa98cYbueXmETRs2JDTTz89oRUdJyfsBXmT6HBmV26edndlO3eZi0B8CRLZyTgQp3bPiLV4vve8zzMRzrusrIxdu3bx+JNP8pclf+Gkk5P5obSEpEaNuHnECMbfOo727dvH7fk74fOZg+Uu4SZYn9eA5Yyqxtxy7rnnajAAqg6LdfqGeKe0tFQ3btyo48aP16YnJ2tG85ba+KQmenJKqv5+3K26ceNGLS0tjfQw653y8nKdkpOr6RkubZvdURufdJIOvnaEvvTlLl26uVCXbi7U/NXrND3DFfD1AdZqFMiIUC2ByprS0lJNS8/Q/NXrKq9dba5fvOG+x9LSM7RtdidNS8/QKTm5Wl5eHumhhYUpObl6dp++lfdB/up12rX3+dqk6cna4pQ2cX/+NRGqZ3Ggcibhp6gSFXdNKu+lplw5sUhKSgr5i5ew+pPPWPD2GhavWc8T73xIdpez+Mtf/sKASy6rMp2QKCQlJTHz3hm8t+JfvPrKUn772zHsL/iWkkMHATOlFwzGqd0Zz0i9eW+viWu/rrKyMvIX5zN5zhNVHPdzHn6Sho0a8cjrK+P6/L1xesb4ymAcLpcQo+AkKD5rUsVBQVFvjOCpjqe/yLBfDmfAwEEkJSXx8949o9p/QkSGiMgWEdkuItX+WSJyg4h8bi//FpHuHtt2isgXIrJeREKTSMvGn09TcdFB0uLwxaEm/EXqxaNfV2FhIWkZLkclN9WVaYfMx+/5e+P0jCm3t3lbWsLlIhESBSdahY7BAEbwOOH0Zv3+2v/QQBpUqa8z9+E5UZP1WUQaAo8DQ4GuwHUi0tWr2TdAf1XtBvwRWOS1faCq9tAQ+wn5cmqfM3kMDRs15owzuySchTDRrFpZWVmUFBc5KrllxUW4WlpZnOP1/KOROis40Sx0AiGep2YMFkbwVMX9Zj3yrgc4eqSMo0eOVFHwADp37hyN01J9gO2qukNVfwReAK72bKCq/1ZVdyLyj4C29TU4z2iasQN7c+vgC+jQuQuLVn6WcBZCSLxUDSkpKYweNZpHp46rouTOnz6RQddcS9Pk5Mp18Xj+nrifp95TUZn1PI5QWHCiWuiAn1Bw4ntqxmBhBE9Vdu3aRYNGjZj+68u5f+wIfjewF08/OJO0zGakZmREs4LXBtjt8b3AXueLUcBbHt8VWC4in4nIGF+dRGSMiKwVkbX79+8PeHDuCtObv9rED6UlPPz3dxkzYxYNk5ISzkIIiZmq4cHZs7j4wvOZMLQfo/r1YNylffl+106uuHEUEP/n78anC4S9vb6MCqGwPTsJnfP8tPcldBR4SlW9rTuAJXSAMQDt27cPaoCe83u+amMY4psHZ89i6rRpTBjaj5OaWuGb6c1bcOv9c4HEETwATy58itbts5n66KIq4eGL77uLfd8XMn/BYzz6yNyomZrywKnKiOPPWUQGYsmafh6rL1TVQhFpCbwjIptVdXW1HVoyaBFYYeLBDrKkpARX8xackt2xyvpETN//4OxZTJueR+6wQdXqmsUjSUlJPDp3LvfNnMm2bduoqKjgmWefY9ovL0uI8w8UBaQejAp1zoMjIr8CBqvqaPv7TUAfVZ3g0HYg8ATQT1UP2uuyPIUOMMFJ6HhSlwJ4dS3THi8kWk4cN2VlZVUEz7PPPht0QclYxl9xwFsHX0C/y39B0d7CWhckDWceHBHpC8xQ1cH299sBVPVPXu26AX8HhqrqVh/7mgGUqarfk6yNrPF3jadcNYCtmzfHZVVtfyRCHhx/JNr5+3vOunPg1OWZG6icCcUUVQHQzuN7W6CajdsWOvnA1W7lBkBVC+2/+7CEUp8QjKlWJJI/zqGSEse8AfGs3IBlNj/nnHPo3bs3C+bPi1qH2nDhz/EzNcPFVb/5XTRPpXwKnC4iHUWkMXAt8KpnAxFpD7wC3OSp3IhIsoikuj8DlwFfhmpgnuUZjMNxdVJSUqLVr6teSPTz9+RQzU1CRigUnKgVOsFi/HESj0QTPP4cP4/98APNT2kTtc7WqloBjAeWAV8BL6nqRhEZKyJj7WZ3A82AJ7wiM1sBa0RkA/AJ8Iaqvl3XMfkqz/DA/fcZh2ODIcKEpFSDiFwOPAo0BJao6v1ugaOqC0UkH7gG+NbuUqGqvUSkE5bVBix/oL+p6v01HS8Qs7GvKZhG/BSL74ln6uhEm64yJBY1lWgo2reX3GGDKNi9K2jFL9FKNdRUnmHPnj2cdnpnZr/8dhWfnLpc41gk0aZoEp1ASjTUxxRV3NaiCsTXxvjjGBIRzxpdjZueTEnxIQZc/StG3XUfJYcOVqufFAyJpOD487VxKy+FhYVcPHgo895eU61/MDW+YhX3vbZ48WLSXJmUFB1i1KhR3P2Hu9i3b59ReBKAcPh71qcPjiFG8U6l3cihdEOi+CQlEu5w5oLdu1j57nJGjxrNR2+/Ru6VA6Iye3G0Ekgiu0TLBeONd0LJh/7+Di8ufYWsNm2rVVw3xCeR9Pc0Ck4C452roILqeQuMT1L8kpKSQteuXXlswfyEc7YOBYEoL4mYC8aNU6mGN55bTPOstix4e03c16YyRJ6EVnB8JgAMU+EvgyFaSTRn61AQqPLimeE4Wmt8hQNvC9fRI0f418vPM3H2/ISoTWWIPAn9mhbvIdEGgyG8BJLIzj0lOPPeGQnlaOtp4bIcsPeQWsOUXjz6IxkH68gRtxYcY50xBIJn/hKDIVg8/ZlqmuJLNCuZt4XL1bI1JYcOJow/kq8UAsbfqP6IWwUnURPZGQLDCB9DKEk05SVQPKfnbr9mMALVasLFqz/StOl5rPjwY6Y98VdmLV1m/I0iQNyGiRtqxjt8LwnL0dibeCzfkJM7lRUffsyvxk8l+8yzOHakrE7h0QaLRAoTNwSOe5qmZcuWzPzjfSxZsiSuS6QUFxfTpl17ECEtsxmlRYe4ePh1XHHjKKb98rKEyX8ULhI+D47B4AsjfMKHUXAMgRDvfikjR43iww1fMvWRp6ol01y3Ynnc5z8KN4HKmfhRmQ2GAMnJzaX9GV2qCZ83nltczdkx3gWxwRAJ3FN68UhZWRkvL13KI6+trBItNmHWPCZdOYCGIlX8jYyMCR9x64NjMDjhFj5u5QZ+Ej7/evl5DhcdIisry/joGAyGWlFYWEi6q5ljtNhJTZvyy1/+DykpKUbG1ANGwTEkFIEKH+8MrMZB0GAwBIK/BJBHy0p5ZO5coHqWZyNjQo9RcAwJRSDCxykDq0lIZggGk34gcfGVAPLx2ycx9ndjycjIMDKmnkgoBce79pKptZR4BCJ8AqkxFAvU9JCNxYewiAwRkS0isl1Eqr3qisV8e/vnItIz0L6hwEw7BI7n/ReL96I/nLJXDzivd2UCyFiXMYE8S6PieeuUKybYBRgCbAG2A3kO2wWYb2//HOgZaF+n5dxzz9XaAKg6LNZlMCQK5eXlOiUnV9MzXNouu5OmZ7h0Sk6ulpeXq6pqaWmppqVnaP7qdbp0c2Hlkr96naZnuLS0tDTCZ+Af9/mlpWdo2+xOmpaeUeX8atpeF4C1GgKZ4rQADYGvgU5AY2AD0NWrzeXAW7bMOR/4ONC+TkuwsmZKTq727Ne/8t7JX71Oe1xwkY6fMLFW1zMe8bz/2nToqCenpGrTk5NDfi9GA6Wlpbply5ZqMiPWZUwgz9JwPm8DlTN1tuCISEPgcWAo0BW4TkS6ejUbCpxuL2OAJ4PoWy8YS46FW+uO58riNWWfjfUCiTXN7cfw3H8fYLuq7lDVH4EXgKu92lwNPGPLwY+ADBE5JcC+dcLXtMP4WfNYtOgpxk+YaCw5VL3/zr14MKf9rEfcFt/0lQAy1mVMrBCKKaqoFjqBYqpmW7grjCdCZXF/2WdjtUBiTXP7e/bsieW5/zbAbo/vBfa6QNoE0rdO+Jt2cLVszXsffhQ3D+7a4nl/NklOSejim7EqY2KJUCg49SJ0RGSMiKwVkbX79++v86ANBn8EU2MomvD3kE1JT+ejjz6K5bl/cVjnnanUV5tA+lo7qKWs8efAfqTkMOPuT4wHtz88789Aim/GM7EqY2rCbe0HyIzwWEKh4NSL0FHVRaraS1V7tWjRIsghGgy1I9ZqDPl7yO7fu4ebbr6ZQwf2c6Dwu2rbY6DgYQHQzuN7W8D7KeirTSB9gdrLGve0w2N5k6pMOyzIm8Sga67llOyOCfHg9kfVCuOtKfVxr8bAvRgyYk3G1EQVi3+ExxIKBadehE4o8FlhPFwHNBjqGV9z+/OnT2TI9bfw6Our6NT1bO79za9jce7/U+B0EekoIo2Ba4FXvdq8Coywo6nOBw6r6vcB9q0zD86eRb9zezDu0vMZd2lfJl81kA5nduXGnDsS7sHthOf9eexIGRcPv4750yfG4r2Y0ATzLK3WJjW1voZZ9ygqrHIPO4CO/BSdcJZXmyuoGtnwSaB9nZbaRlF5eGCbaCof4HktzDWKSSqjVDIyNKN5S01OS9dhI8fqS1/uqozUaJqcrGnpGY5RZHWBMEZR6U9RUluxIqLutNeNBcbanwUrcOFr4Augl7++NS21lTW3jp+gXc/trY+9/UHlNe/Zr79Oycmt1f7iiaKiIv3NyJGalpGhbTt01OTUNG2abEVRhfJeNNQf9f28CFTO1HmyT1UrRGQ8sAwrFHOJqm4UkbH29oXAm7Zw2Q78APzGX9+6jqkmXKmpiIOzbL1qllGEZ1XxJKwnhPuvN4l6jWIJ99z+TTfewNCrhvHYsg9ompxcud3VshXNW7Tin6+8THJyckzVwFHVN7Hkiee6hR6fFbg10L7h4tFH5jJteh53XntltarZiUpFRQXTpuexePFi0lyZ6AnlkkEDeGTuXJKSkkw9JkPIMdXEDYiIo+OTALF4fxgsysrKaNO2HXNfe6+KI2fRvr3kDhsUlqrpppp4VUwhxZ/IyZ3Kqk/WVkbwuaej+vfpxdyH5/jsZ65h9FPfz5BA5UxCZTI2GBIJk2sj8sSbA2ltqU1pApMVOnbw6ZMTYYu/UXAMQRMVKbgNAWFybRiigdqUJojhhJQJx6GSEkcfmEMlJREdl5miMgRtXjRTWrFHfZn5zRSVwYlgp0sjMb1qiB3MFJXBYKjETJUYIkmw06WxXozSEB3EdspEQ0gwUWUGgyHcPDh7FtOm55E7bFCNkWVVEwJWteAkei4hQ+AYC44haudPDcFTVlbG1q1bE7ocgCE6CaY0gdvis2D6BL767BOOHjliHOQjTCz6XhoFx2CIA0zEiSFWCGS6tKKighN6gq0b1jFn0hhGXvAzJl3Rn5/36mkc5COEuxBzLBVgNgqOISjc2rp3OGAjwJWSYqwHEcJEnBjiiWnT83j/0/+w4O01LF6znife+ZDO3XrQoEEDDhw4wLJly9izZ0+kh2mIckwUlSEo/EVQpaVnkOZyUXzoILfcfAuPzH045ivjxgLBRJyEO5rKRFEZ6oqv+/lA4XdMunIAx49XkJKWQVlJMd279+D9VStp0qRJBEecGNQUPeuZEd8TV2pqyN0dTBSVod6xrAcf8Ojrq3j3gw/pe2E/M0VSDwQScWKmsAyxgq/7+Y3nFpN9ZleefOcj8t9fx5PvfMSRE/Dz/gMiM1BDFaJxCssoOIaQ4ZmhdOqji/jii8+ZPCUnwqOKf6yIk4OV4bduivbt5XDxIbKysswUliFm8IygcnP0yBHe/b+/MfXRRVXkTO7chXy+YYOZrjI4YhQcQ1hwtWxFRvOWPP3008Ynpx5o374DD036bZUcIw/njKVz5zMAgk6TbzBECqecOTs3b6RJ05MdrZTJaels2LAhEkNNKNwFmL2XaHZCiOaxGWKYon17OVJymHRXJoWFhXTu3DnSQ4o73NWZ8xfn07DRSRw9UsbEyy8iNbMZRw4Xc+HQYXz41qts27atxiks8/8xRBPeOXOKiw5y7Nixanlxvt/5DWXFRZx22mkRHG1849O3BjiEpeREK8aCYwgKX0XVGjdoUMV6sCBvEuddMpTDRYdIi+I8CbGMe9rpkddWsnjNeh5f9gGndTuHs/v0ZeGKTxkzYxbprkyAaiZ/iK2kaSKSKSLviMg2+6/LoU07EXlPRL4SkY0iMslj2wwR+U5E1tvL5fV7BoZg8M6ZU1hQQPfuPXg4ZyxF+/ZyvKKCRTPyyPnFxaS6Mul5bi/jUxYmfPrWRHRUgVEnBccIncTDKSlgeXk53Xr15tbBFzDu0r5MunIAZYeLWfPGP2iaksoZZ3YxwifE+KrOPHHWPD5+5y3gJwXm9NNPj4eq4nnAv1T1dOBf9ndvKoBcVe0CnA/cKiJdPbY/oqo97OXN8A/ZUFc8c+a8v2olyQ1g3KV9GXlhN77d+hVPLP+QP6/+j/EpiwKisaJ4XS04RugYSEpK4sMP1jB61GiOFBfRqFEjGp10Eo8v/zcL3/vUCJ8w4C9yKiXDxc7NG6soMHFQVfxq4Gn789PAL7wbqOr3qvof+3Mp8BXQpt5GaAgrTZo04dOPP2LTxi/RE8erORwbn7L6x1OBicaM+HVVcIzQMQCWkvPYgvls27qFivIfjfAJM06RJmBZZg7t/Z4Hx91SRYEJJk1+lNJKVb8HS6YALf01FpFs4BzgY4/V40XkcxFZ4mRt9ug7RkTWisja/fv3133khpBy/PhxXM1amI7YbmIAABHhSURBVEKcUUCkFZiaqKuCY4SOoQolJSVkZDY3wqeO7Nmzx2+2Vl/VmR+/fRI3XH893xXsdlRgormquIi8KyJfOixXB7mfFGApMFlV3dL3SeBUoAfwPfCwr/6qukhVe6lqrxYtWtTybAzhwp9yHys+ZZGksUM9KRGhsUSzu3DtqFHBMULHEAw1CZ+0tDS2bt3Knj17TFkHB44dO0bv886nfXY21910M+2zs+l93vkcO3asWlunaacB5/Vm0VNPRaUCUxOqeomqnu2w/BPYKyKnANh/9zntQ0QaYcmZ/1XVVzz2vVdVj6vqCeDPQJ/wn5EhHPhS7j2nZGOxMGR9UU51h2G11zsRjb41gVKnUg0isgUYoKrf20Jnpaqe4dCuEfA6sExV5/rYVzbwuqqeXdNxTfr06CYndyqrPllb6QDrFj4NfjzGls2badi4MUdKLUvPf38oY9So0Tw4e1YsTZeEjd7nnc+RE5A7d2HltXs4ZyzJDeDTjz9y7BPu8gvBEK5SDSLyEHBQVWeJSB6QqarTvNoI1lT5IVWd7LXtFLe1WUSmAOep6rU1HdfImujEnSJhyZIlpGW4KCkuovzIEX4o/7Gyjb+yAolMTSUXYoFA5UxdFRwjdAzVcBI+p3fuzPFGJ9GyXTZ7du1kwqx5VZSf/n16MffhOZEeekTZs2cP7bOzefKdj6rVlBp3aV++3fkNrVu3juAIayaMCk4z4CWgPbAL+JWqHhKRLCBfVS8XkX7A+8AXwAm76x2q+qaIPItlKVZgJ/A7t+zxh5E10Y2ncp+amlr54BaMguMLo+AEfhAjdAw+cQuftLQ0zjizCw+8+AbTf305815fWWNRyERk2bJlXHfTzeS/v67attE/P4fnn32awYMHR2BkgWOKbRoiheeD2yg4vkkkBadOcwKqehC42GF9IXC5/XkNPpIdqupNdTm+IbpxO7Ru3bqVNFcmqidINRl1fdK9e3fKSoqrZWt1Z4Xu3r17BEdnMBgMsYXJZGwIO27HY5EGlJroB5+0bt26SrZW+KmmVLfu3aN+espgMEQ/jXCuKdUokoMKE0bBMYQdd9TDkvvuoN/lV7Mgb1IsZ9QNK57ZWkf//Bx+f8n5NPzvUd564/VID81giBlcOD/EYyHyJ9z86JGEz/N6lEPcRZsZBcdQL7hDmj9861V2bf2KcZf2ZezA3uQMGxhrGXXrRFlZmd/weHe21q+3b2NQ/4to0qQJxaVlnHra6abchcHgB89w5iKv9dGQVTfc1CY03medKYfimrGIUXAM9YI7k+53Bbv56IM1fLvzG9asXMF3u6snpPOnBNSkIEQrFRUV5OROpU3bdlw8eCht2rbzq7DMfeRRvtmzj0deX8m8t9eYchcGQw0EUiognvPjxLuyUhuMgmOoV9yOx61bt66WUdefEhCsghBtuCt/z33tvRoVFl+FNE25C4OhbtSkBMSzApSImMxqhqjBUwnwzJHjVgJ8bYv2/DluhcU9dvhJYckdNoiZ986oouj5K6Rpos0MhvDhVoC8kQS2gsQyxoJjiBie003+rBaLFy8mf3F+1Fg0gp0mC0Rh8dxvWlqaqbVjMNQRJ2sMQGaEx+WPcFmQfO033i0cRsEx1DtO000TJ00iNcPlqAScnJJCmo9t9VnAs7bTZDXV52rZsmWV/Z5xZhc6n3EGj99uos0MhtriczoqoqPyT7j8aHztt4L4jjYzCo6h3nHyR1m/ZTtFBw84KgE/lJVRUlwUcYtGMH40ntRUHHDmH++rtt8TjZvQsPy/VQppJlK0mcFgCI7aFsWsyTE7lqlTqYZIYdKnxy5lZWW0aduuij8KWA/8CUP70bnbOUyYvaBanSrAsYBnfdWw8jfuQMpMONXnGjlyJHf/4S46ZHf0ud/NX22ipKQkKgppBoIp1WCIFvyVJPDGlZrKoZKSiJcxCNfxI31eoaZeSjUYDMHizx+lWYtWnHPm6eQOG0RahovDRYf45S//h7v/cBcpKSlMm55Xuc2tILgtGuGuqF1Xx193mPzMe2dUGae7jIWv/ZaUlBiHYoMhxPh6qLtSUx0dipOg0ofHs228WDriFTNFZahX/PmjlBYXsWD+fHZ+s4NLBg1AUd59bxUdsjsybXoeD86eRcHuXby77C0Kdu+qtNx4+8WMnzCRTZs2hdT52HvcR48cofCbr/l+5zdBTZO5w+TdSlhN/jnGodhgCA9OjrdFpaVVEgO6FaFyLJ8Vl0f/otLSsISRZ1J9msk9XkNwGAXHUK/U5I+SkpLCzD/ex4atX/PIa9WT3HkrCJ5+MXNfX0nfocPIz89nwCWXhTRXjnvcj98+iUUz8vjdwF7cP+ZGcn9xCad37kyTJk3Cdj0MICKZIvKOiGyz/7p8tNspIl+IyHoRWRtsf0P8UJNPSm0ceosc2ocymZ4rNTUsx6itf07M4+RgFOiCpWy+A2yz/7p8tNsJfAGsB9YG2997Offcc9UQu5SXl+uUnFxNz3Bpu+xOmp7h0ik5uVpeXq6lpaWalp6h+avX6dLNhZVL/up1mp7h0tLS0sr9eLcdNnKsdr/gosrv+avXac9+/XVKTm6NYyotLdUtW7ZU2b/TuHv1OU/P7Nm7VseozfWINTx/36FcgAeBPPtzHjDbR7udQPPa9vdejKyJXwBVh8V6LDq3C7SPJ67UVCd9RV2pqXUaVyITqJypk5OxiDwIHFLVWSKSZyso0x3a7QR6qeqB2vT3xjj+xQdOfjNbt27l4sFDmff2mmrtJw/px7vL3qr0SfFse/TIEX43sBfzXl8ZlBNwcXExU3JyWLr0FdIzm1FSdIhRo0bx4OxZVcpHuMdbk6MxUGtfoHD7EdUH4XIyFpEtwABV/V5ETgFWquoZDu124ixrAurvjZE18Uugjree7QSCctbNTEtztLy4sK1BDn38jSsJK7S7yr4S0BcoUDlT1ymqq4Gn7c9PA7+o5/6GGMZ7ugmC80nxbFu0bw+pASTTc+POaZPVpi2vvvkWx1XpOegyHvr7Oz5Dv/05GjduejK3jp9Qp1ISTtfDUEkrVf0ewP7b0kc7BZaLyGciMqYW/RGRMSKyVkTW7t+/P0TDNyQioc7FU+G0L5Nl2Sd1VXCM0DGElGB8UjzbijSgNADFyJ0teEpOLis//pQFb6/hz6v+w/w3VvHt5k288dxinxmS/SpfRYdYvuK9ylw2D7z4BstWrWbylJxwXKa4RETeFZEvHZarg9jNharaExgK3CoiFwU7DlVdpKq9VLVXixYtgu1uiDM8/Vf84StzcqjxdkIGTL0sX9Q0hwW8C3zpsFwNFHu1LfKxjyz7b0tgA3CR/T2g/t6LmRePb4LxSfFs62rewqd/jLtdWnqGZrXP1kYnneTo55OSnqHPfbZN22V30i1btlQ73vgJE6sdo/sFF+nga0doclq6PvPJVzps5FhNTkvX1u06aKOTTtLxEybGpD9NbSB8PjhbgFPsz6cAWwLoMwOYWtv+amRNXBOsb0xNffDlN+NjPT58anwew7NfgvvoBCpnjNAxRC2BOP56tt24caOOnzDRUTGakpOrPfv11/zV63TBW+9r63Ydqig37qV1h4563//+o5pDs5uNGzdqRvMWmpKeoa07dNSU9AwdNnKsvvTlLm3doaMOuubaao7OPS64qE5OyLFEGBWch6jqJPygQ5tkINXj87+BIYH2d1qMrDEESqgUHF/7UqPgVBKonKmrk/FDwEH9yUk4U1WnebVJBhqoaqn9+R1gpqq+HUh/J4zjn8Ef3s663s7B/hySJ105gNPO+hmD+p7nmCHZva8HXnwD1RO4WramaXKy1feK/igw/41Vtcp2HA+E0cm4GfAS0B7YBfxKVQ+JSBaQr6qXi0gn4O92lyTgb6p6v7/+NR3XyBpDoDg5B/tzSvbnHFxTFuZgHJ3jkfpyMv7/7d1bqFzVHcfx74+oIMdLtCHxRMULhGJ9kIYSjYpIY0TOS+oNLBS1lKhIqlX6cFDxyYcoxN5sH8QIEaR9UWvQeEcpItqL5CRKaqNBMZyQxMQafTGCfx9mHzk5mcueM3v2rFnz+8Aw+8zs2fu310z+rOxZs9d6YLWkncDq4m8kLZW0pVhnCfCmpCngn8DzEfFiu9eb9WLuYN25g4OPHxtj1XU/50+TR05mueGuWyGCn668sOWcTzPjfh5/4B6OHzvh+87N7397O2MnnZzEpKA5iogDEbEqIpYV9weLx6cjYqJY3hURFxS382c6N+1eb9ZPp9B6Mst2v3xqdd0aTz3QnZ7aKyIOAKuaPD4NfF90gAu6eb1ZlWYPDp7pfPzi7nvY+MB93L56JYvHxzn0/8+59ppreOvVl1i4cGHb7T304Pqjpo24+Zc3c/jrwzy28bEj9gO+KrHZqDrI/M6stOv89Gvwco7cIbTszf611cxknYcOHmD/7k9Yu3Ytd/x6XVfXnmk1rxTAggULeGTyTtat/8MRk4L6qsRmeWs1j1XVVwuuaz85cAfHRkKzsy4zk3XOvaBfWTNfhc32u4c3tJ0U1MzyVNfF9kbton696GmQ8aB44J/NV11XC87hqsTz0a9BxoPiWmOWnrJ1xmdwbKQ0O+syzPsxM7PmPJu4mZmZZccdHDMzM8uOOzhmZmaWnaEcZCxpP/BJD5tYBHxWUZxhzgBp5EghA6SRI4UMMP8cZ0VENjNU9lBrhv19zC0DpJEjhQyQRo5eMpSqM0PZwemVpH8P+pceKWRIJUcKGVLJkUKGlHIMq1TaL4UcKWRIJUcKGVLJUUcGf0VlZmZm2XEHx8zMzLIzqh2cRwcdgDQyQBo5UsgAaeRIIQOkk2NYpdJ+KeRIIQOkkSOFDJBGjr5nGMkxOGZmZpa3UT2DY2ZmZhlzB8fMzMyyk30HR9L1kt6X9K2klj9Jk3SVpA8kfShpsg85TpX0iqSdxf0pLdb7WNJ2SVslVTLLX6djU8Mfi+e3SVpexX7nkeNySV8Ux75V0v19yPC4pH2S3mvxfF1t0SlHHW1xpqTXJe0o/o3c2WSdWtpj2LnOpFFnXGO6ypF/jYmIrG/AecAPgTeAn7RYZwHwEXAucBwwBfyo4hwPAZPF8iTwYIv1PgYWVbjfjscGTAAvAAIuAt7pw/tQJsflwHN9/jxcBiwH3mvxfN/bomSOOtpiHFheLJ8I/G8Qn40cbq4zg68zrjFd58i+xmR/BicidkTEBx1WWwF8GBG7IuIw8DdgTcVR1gCbiuVNwM8q3n4rZY5tDfBENLwNLJQ0PoAcfRcR/wAOtlmljrYok6PvImJPRLxbLH8J7ABOn7NaLe0x7FxnkqgzrjHd5ei7QdeY7Ds4JZ0OfDrr790c/Sb0aklE7IHGmw4sbrFeAC9L+o+kWyrYb5ljq+P4y+5jpaQpSS9IOr/iDGXU0RZl1dYWks4Gfgy8M+eplNpj2LnO9Pf4XWO6l3WNOaaKjQyapFeB05o8dW9EPFtmE00e6/r38+1ydLGZSyJiWtJi4BVJ/y164vNV5tgqOf4KcrxLY46RryRNAH8HllWco5M62qKM2tpC0gnAU8BvIuLQ3KebvGQkry3hOtM+VpPH6q4zrjHdyb7GZNHBiYgretzEbuDMWX+fAUxXmUPSXknjEbGnOP22r8U2pov7fZKeoXHatZfCU+bYKjn+XnPM/uBHxBZJf5G0KCLqnBSujrboqK62kHQsjcLzZEQ83WSVJNojBa4zbaVQZ1xjujAKNcZfUTX8C1gm6RxJxwE3AJsr3sdm4KZi+SbgqP/xSRqTdOLMMnAl0HQEfBfKHNtm4MZiNPtFwBczp7kr1DGHpNMkqVheQePzeaDiHJ3U0RYd1dEWxfY3Ajsi4uEWqyXRHplwnenvZ8k1pgsjUWOqGq2c6g24mkYP8WtgL/BS8fhSYMus9SZojPD+iMYp56pz/AB4DdhZ3J86NweN0f9Txe39qnI0OzbgNuC2YlnAn4vnt9PiVyA15FhXHPcU8DZwcR8y/BXYA3xTfC5+NaC26JSjjra4lMap4G3A1uI2MYj2GPab60wadcY1pqsc2dcYT9VgZmZm2fFXVGZmZpYdd3DMzMwsO+7gmJmZWXbcwTEzM7PsuINjZmZm2XEHx8zMzLLjDo6ZmZll5zuIOlk+H+nDggAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 576x216 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 3))\n",
    "\n",
    "\n",
    "ax1.scatter(X[y_km == 0, 0], X[y_km == 0, 1],\n",
    "            edgecolor='black',\n",
    "            c='lightblue', marker='o', s=40, label='cluster 1')\n",
    "ax1.scatter(X[y_km == 1, 0], X[y_km == 1, 1],\n",
    "            edgecolor='black',\n",
    "            c='red', marker='s', s=40, label='cluster 2')\n",
    "ax1.set_title('K Means Clustering')\n",
    "\n",
    "\n",
    "ax2.scatter(X[y_ac == 0, 0], X[y_ac == 0, 1], c='lightblue',\n",
    "            edgecolor='black',\n",
    "            marker='o', s=40, label='cluster 1')\n",
    "ax2.scatter(X[y_ac == 1, 0], X[y_ac == 1, 1], c='red',\n",
    "            edgecolor='black',\n",
    "            marker='s', s=40, label='cluster 2')\n",
    "ax2.set_title('Agglomerative Clustering')\n",
    "\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Clustering using DBSCAN\n",
    "\n",
    "Unlike K Means or Agglomerative Clustering, DBSCAN is a density-based approach to spatial clustering. It views clusters as areas of high density separated by areas of low density. This approach has several advantages; whereas K Means focuses on finding centroids and assoicating data points with that centroid in a spherical manner, the DBSCAN algorithm can identify clusters of any convex shape. Additionally, DBSCAN is robust to areas of low density. In the above visualization, we see that Agglomerative Clustering ignores the low density space space between the interleaving circles and instead focuses on finding a clustering hierarchy that minimizes the Euclidean distance. While minimizing Euclidean distance is important for some clustering problems, it is visually apparent to a human that following the density trail of points results in the ideal clustering. \n",
    "\n",
    "For more information on the DBSCAN algorithm and its implementation in scikit-learn, check out this resource: http://scikit-learn.org/stable/modules/clustering.html#dbscan"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cluster import DBSCAN\n",
    "\n",
    "db = DBSCAN(eps=0.2, min_samples=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, let's fit our model to the data and generate predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_db = db.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, let's visualize the model applied to our data. We see that the DBSCAN algorithm correctly identifies which half-circle each data point is generated from."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XuUVOWV9/Hv1sbBobulmzs0iCRiRAeQtHhLIsEYlJlAZkbf5SVecQghKggE0KylyCQuUCBqvDAOkDFmJsbXZN6QRMWoUVecJIqxvSDXMAoVjCC00m000LrfP7qaVHfXvU5Vnar6fdaqRVfVqXOeorprn+d59tmPuTsiIiJhc1ixGyAiIhKPApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiIRSVbEbkEzfvn19+PDhxW6GiIgE6MUXX3zH3ful2i7UAWr48OGsX7++2M0QEZEAmdmb6WynIT4REQklBSgREQklBSgREQmlUM9BiYgU28GDB4lEInz44YfFbkrJ6dmzJw0NDfTo0SOr1ytAiYgkEYlEqKmpYfjw4ZhZsZtTMtydvXv3EolEOOaYY7Lah4b4RESS+PDDD+nTp4+CU4bMjD59+uTU81SAEhFJQcEpO7n+vylASUVpbW1ly5YttLa2FrspIpKCApSEXhBBpa2tjTlz5zGkYShnTTqXIQ1DmTN3Hm1tbQG2VKRdIU6EFi1axLJlyzJ+3bvvvss999yT8/HvuusuPvnJT2JmvPPOOznvLx4FKAmtIIPK/AULeeb59az42a+447Ffs+Jnv+KZ59czf8HCPLRcKlUpnAhlE6DcnY8//rjTY2eccQZPPPEERx99dJDN60QBSkIrqKDS2trK6tWrmXnL7dT1HwBAXf8BzLzldtasWdPpLFdDgJKLfJ0Iff/732f06NGMGTOGSy65pNvzEyZMOFQW7p133qGjhumGDRsYP348Y8eOZfTo0WzdupWFCxfyhz/8gbFjx/KNb3wDgNtuu42TTz6Z0aNHc9NNNwHwxhtvcPzxxzNz5kzGjRvHzp07Ox3zpJNOIt+1UhWgJJQyCSqp7Nq1i9q6+kP76VDXfwC1vevYtWtXSZz5SrgF+Tsba8OGDXz729/mqaee4uWXX+aOO+5I+7UrV65k1qxZNDU1sX79ehoaGliyZAmf+MQnaGpq4rbbbuPxxx9n69atPP/88zQ1NfHiiy/y7LPPArB582YuvfRSXnrppbz2lBJRgJJQSieopGvw4MHsb95H8+63Oz3evPtt9r/bzODBgzUEKDkL8nc21lNPPcV5551H3759Aaivr0/7taeddhq33HILS5cu5c033+TII4/sts3jjz/O448/zkknncS4cePYtGkTW7duBeDoo4/m1FNPzardQVCAklBKJ6ikq7q6mmnTpnHPDbMP7a9599vcc8NsrrzySoC8nPlKZQnydzaWu6dM166qqjo0RxR73dFFF13E2rVrOfLII5k0aRJPPfVU3P1ff/31NDU10dTUxLZt25g2bRoAvXr1yqrNQVGAklBKFVSqq6sz2t+tS5dw5vhG5k6ZyOxzPsPcKRM5c3wjty5dkrcz3yBoTqx0BP072+Gss87ioYceYu/evQDs27ev2zbDhw/nxRdfBODhhx8+9Pj27dsZMWIE1157LVOmTOGVV16hpqaGlpaWQ9tMmjSp04nYH//4R3bv3p1VW4MWSIAyszVmttvMXkvwvJnZnWa2zcxeMbNxQRxXSk+6X7itra1cNe1KTj9pTNygkqmqqipWLF9GZOcOnlj3KJGdO1ixfBlVVVV5O/PNhebESlOyE6FsnXDCCXzzm9/kzDPPZMyYMcyZM6fbNvPmzePee+/l9NNP75Ty/aMf/YgTTzyRsWPHsmnTJi699FL69OnDGWecwYknnsg3vvENvvjFL3LRRRdx2mmn8Xd/93ecd955nQJYInfeeScNDQ1EIhFGjx7NVVddlfV7TMjdc74BnwPGAa8leH4y8ChgwKnA79LZ76c//WmX8nDw4EG/bs5crz2qtzcMH+G1R/X26+bM9YMHD6bc7utXX+MbNmzwlpaWvLXvujlzfdxnzvRVz77kP960y1c9+5KP+8yZft2cuSlf29LS4ps3bw60fbm0R4L1+uuvZ/yafPxOlKp4/3/Aek8ntqSzUVo7guFJAtS/ARfG3N8MDEq1TwWo8pHuF26xvpg7AuNRvet86PARflTvurgBNN5rUgXdTLW0tHjtUb0P/R903FY9+5If1btOX3oFlk2Akr/KJUAVag5qCBCbRB+JPtaNmU03s/Vmtn7Pnj0FaZzEF9T8R7rpt/lI0033PSQbAkwkX5l/YZ4TEymkQgWoeCkoHm9Dd7/P3RvdvbFfv355bpbEE/T8R7pfuEF+MWf7Hqqrqxk5cmTKCe18XfMC+csGEyk1hQpQEWBozP0GQKeBIRV0zyDdL9wgv5jzfV1TPns5+coGEyk1hQpQa4FLo9l8pwLvuftbBTq2pKFjKOxPf/pT4D2DdL9wg/pizmfvpkO+ezn5yAYTKTVBpZn/EPgNcJyZRcxsmpnNMLMZ0U0eAbYD24B/B2YGcVzJXdehsGOPO47DevSgtr5Pp+1y7Rmk+4UbxBdzPno3Xeey8t3LyWZOTKTspJNJUaybsvjyL17W3KfGneyTLrg0Lxlk6abf5pKmG2QWXLJMvWwy/6T0hDWL76abbvLbbrst49c1Nzf73XffnfPxL7roIh85cqSfcMIJfsUVV/iBAwfiblcKWXwSQomGwubdfh9P//T/8tYb/wsE2zNINwkh3e0SvTao3k2yuax4vZzFNy9i+/btnTITVQmictTX1mJm3W71tbXFbtohQS23cfHFF7Np0yZeffVVPvjgA1atWhVkMwGVOqpoO3bs4G/+tleCobB6Fp5/bsnOfwQxVJjuXFZ1dTUjRozgxpsWxWQNNnDyKacyeEiDKkFUkOaWFhy63ZrTqMyQTBiX25g8efKhADx+/HgikUhO7zEeDWhXsHvuXXlooj82SDXvfpsDH/yZrVs2s3//fgYPHlxymWMdvZvFNy9i165dWb2HdOayRo4cCXTuadX1H0Dz7rdZPmcGp0+eyvRFSw714OYvWMiK5ZmvgpqN1tbWrN+7hEfHchvPPfccffv2jVuLL5GO5TYuvvhiDhw4wEcffcSSJUt47bXXaGpqAui03Ia7M2XKFJ599lmGDRvG5s2b+d73vpe0x3Xw4EEeeOCBjJYBSZd6UBWqtbWVBx54gAlfPp/vLpzVaShs2ezpXHLJJQwcODDrYbawyGWoMN1MvUQ9rbkrVvLco2v54P33C1odXXX8ykvYl9uYOXMmn/vc5/jsZz+bwbtKjwJUheroHVx03fX0HTyEWf8wga9POoPZX/o8b+94g6/N+Gqxm1h06c5lJetpVfeuo3n3nw7dL0QlCK1tVV7cw7vcxs0338yePXtYsWJFpm8rLQpQFap///7se2cPX/vCKWx44be4O6MaT+Hb//VTPm5rY9iwYcVuYiikM5eVrKfV+m4zdf0HHrqf70oQhbgGTAorrMttrFq1inXr1vHDH/6Qww7LTyhRgKpQi//1W4wYdSJ3/uIZ7l73HHf+4hne2fVHls68XNUKYqRzPVKintbyOTM449wpHNmrV8EqQaiOX3HV1dRg0O1WV1OT9T7DutzGjBkzePvttznttNMYO3Ysixcvzvo9JmLtKenh1NjY6B2ZKRKc1tZWhjQMPTSh36F599tcc+5n2BWJ0Lt37yK2sPS0tbUxf8FC1qxZQ23vOt57dx8jRx7Hls2bOaqunv3vNnPllVdy69Ileb3YNtlnO3fKRCI7d+jkI0MbN27k+OOPL3YzSla8/z8ze9HdG1O9Vll8FSjZWXbffgPYvXu3AlSGEmUNFjqTLrY31zHMpzp+UqoUoCpQ7JxJ17PsQlfLLrdU6I6swUT3C+HWpUuYv2Ahc6dMpLZ3Xafem0gp0RxUBQpDtWylQueP6vgFL8xTIWGW6/+bAlSFKna1bKVC518u14DJX/Xs2ZO9e/cqSGXI3dm7dy89e/bMeh9KkqgAyYbRijHEpol8KSUHDx4kEol0ur5I0tOzZ08aGhro0aNHp8eVJCGHMstWr15NbV09+5v3MW3atE6ZZMWYI8mkhJBIsfXo0YNjjjmm2M2oSBriK2NhHUbTkuYikg4FqDIV5ooCYUjSkOS0TIiEgQJUmQp7RYFiJ2lIfMqulDBRgCpTYR9GUyp0OIV1WFgqUyAByszOMbPNZrbNzLr9JpvZUWb2MzN72cw2mNkVQRy30iUbhimVYTSlQodHmIeFpTLlHKDM7HDgbuBcYBRwoZmN6rLZ14HX3X0MMAFYbmZH5HrsSpXuMIyG0SQTYR8WlsoTRA9qPLDN3be7+wHgQWBql20cqLH2RU2qgX2ABrWzlO4wjIbRJBNhHxaWyhNEgBoCxC5WH4k+Fusu4HhgF/AqMMvdPw7g2BUnm2EYDaNJOkplWFgqRxABKt5Sj13LU0wCmoDBwFjgLjOrjbszs+lmtt7M1u/ZsyeA5pUXDcNUrkKkfmtYWMIkiAAVAYbG3G+gvacU6wrgJ95uG/C/wKfi7czd73P3Rndv7NevXwDNKy8ahqk8hUz91rCwhEkQAeoF4FgzOyaa+HABsLbLNjuAswDMbABwHLA9gGOXvEzPijUMU3mKkfqtYWEJg0CKxZrZZOB24HBgjbt/28xmALj7SjMbDPwHMIj2IcEl7v6DVPst52Kx6dTJS/XajtVbC7VaqxReroV1y229LSkP6RaLVTXzIpkzdx7PPL++26qnZ45vZMXyZWntQ18+5W/Lli2cNelc7njs192em33OZ3hi3aNxC+vmcgIkkm/pBihVkiiCoC6I1DBM+ct2zlEVIaQcKEAVQUcmXs9e1ez63z/wwfvvA8rEk+6ymXNURQgpF+rrF0H//v3Z984epk/4NLX1fWhp3sdZ513I339lmjLxpJtbly5h/oKFzJ0ysducYzxab0vKhQJUESz+128xYtSJzF52z6H5pzsXXMvNV/wfZeJJNx2p34tvXpTWnGPssGDXxAqdAEkpUZJEgSXLyrrm3M+wKxKhd+/eRWyhlINMknCUbCOFpiSJkEo2/NK33wB2796d0f60sJzEk05FCK39JGGnAFVgQVWC0JeLJJNORQhl+knYKUAVWFCVIPTlIulIdCmCMv2kFChAFUGuBTn15SK5UtFhKQUKUEWQa0FOfblIrlR0WEqBAlQRZVsJQl8ukisVHZZSoABVgvTlIkHQ2k8SdroOqkSporkERddBSaGpmnmF0JeLiJQaXahbIVTRXHKli70lrBSgRCpUqou9Fbik2DRZIVKhYi/2jq3XN2/+fA6zw7TYoRSd5qBEKlCqosUjR5/ENUu/m/VqzyLJaA5KRBJKdrH33xzZi/OvnqcqJVJ0gQQoMzvHzDab2TYzi1sMzswmmFmTmW0ws2eCOK6IZCfZxd5/bt3P8E+d0OlxVSmRYsg5QJnZ4cDdwLnAKOBCMxvVZZvewD3AFHc/ATg/1+OKSPYSXex99/WzqKrqwYfvd+4pqUqJFEMQPajxwDZ33+7uB4AHgaldtrkI+Im77wBw98wWPSphyoSSsIpXSWLCKSfzL1ddpSolEgpBpOQMAXbG3I8Ap3TZZiTQw8yeBmqAO9z9+/F2ZmbTgekAw4YNC6B5xdFR6UGZUBJWiZaS7/jdnTtlYrcqJSKFlHMWn5mdD0xy96ui9y8Bxrv7NTHb3AU0AmcBRwK/Af7e3bck23epZfHFVnW48aZFaS+5LRJGqlIi+ZJuFl8Qp/IRYGjM/Qag60xqBHjH3d8H3jezZ4ExQNIAVSq69pbea97LRx99zHfWPtUtE2rulIksvnmR/uAl9DqqlGRKgU2CEsQc1AvAsWZ2jJkdAVwArO2yzU+Bz5pZlZn9Le1DgBsDOHYodF3d9js/e5phxx3PL36wutN2yoSScpaqMoVIpnIOUO7eBlwNrKM96Dzk7hvMbIaZzYhusxF4DHgFeB5Y5e6v5XrsMEi0uu287/wbTz78Qz54//1D2yoTSspZ1xO1FT/7Fc88v575C+JeeSKSUiDXQbn7I+4+0t0/4e7fjj620t1Xxmxzm7uPcvcT3f32II4bBskveDySNzZtAJQJJeUt0YmaLvCVXCidLEexFzx2LRnzQWsLt868nKPq6pUJJWUt2Ylax7B2NvNZUtlU6ihHyVa3nfHVGfwxspMn1j1KZOcOVixfphRzKUvJKlNoWFuypQAVgGRLZ2u9JqkEyU7UNKwt2VI18wApvVYqWcflFmvWrOl2ga9GDiSWlnwXkaKIPVEDdNIm3Wi5DREpiurqakaMGMGNNy3SNVGSEwUoEQmcromSIGiIT0QClWy13rlTJhLZuUPDfRVOQ3wiUhTpXBMlkg4FKBEJlK6JkqAoQIlIRlItwqlroiQoClBSkuprazGzbrf62tpiN61sZVKtPNnF6yLpUpKElCQzI95vrgFh/p0uZXPmzst4EU5dvC7x6ELdgOgPLJxiA1Q90Bxnm7qaGvbt31/AVpUvZeZJkJTFlyMtvlY6mgGPc2tuaSlms8qKMvOkGBSgEkh0oeF1c+Z2myBONWksUuqUmSfFoCG+OJINZ8w8+1T6DhhE63vvcsUVV4DB99Z8j9q6evY372PatGkqjlkAsUN8BpqPKoBs5qBE4kl3iE/fonEkG86o6z+Q6//tAY7sVc2Nl/wT9QMGHgpkHX+w8xcs1B9sntXV1GAawiuoW5cuYf6ChcydMrFbtXKRfFAPKo5kPajZX/o8K596AYDpEz7Nnb94RpPGBVRfWxt3bqkO2NflMfWg8kOJQ5KrgiZJmNk5ZrbZzLaZWcJqkGZ2spl9ZGbnBXHcfEl0oeF3F85i4j9fwJG9etG8+0/U1vfRpHGBNbe0xE+IoD0gxd7qamqK1cyypkU4pVByDlBmdjhwN3AuMAq40MxGJdhuKbAu12MWQuyFhrPOOYOZZ5/KwGHD+cqcGwCo6z+Q/fv2atI4RNy9000p5iKlLYge1Hhgm7tvd/cDwIPA1DjbXQP8GNgdwDHzrqqqihXLlxHZuYMn1z3Gv/zLdPZE3mT/vr0AfPh+K73r+3D7vJkq5yIikgdBJEkMAXbG3I8Ap8RuYGZDgH8EJgInJ9uZmU0HpgMMGzYsgOblprq6msGDBzPzazPg3pWdJogvv/xyMDRpLCKSB0EEKIvzWNeZ6duBBe7+kVm8zWNe6H4fcB+0J0kE0L6stbW1MX/BQlavXn0ojfzSSy/lazO+yrBhww71kr61eLEmjaWiKXFC8iGIIb4IMDTmfgPQNUOgEXjQzN4AzgPuMbMvB3DsvIp3se7/vPQyq1av6fRHqEnj/IotDFtF92QIJUQUjyquSD7lnGZuZlXAFuAs4I/AC8BF7r4hwfb/Afzc3R9Ote9i1uJT7bHwUGHY8NLFu5KNgqWZu3sbcDXt2XkbgYfcfYOZzTCzGbnuv1hUe0wkudbWVlavXn0oOEH738fMW25nzZo1Kv0lOQukkoS7PwI80uWxlQm2vTyIY+ZbbO2xrj0opZGLpHcSN3LkyCK1TsqBisUmoFVBRZJTAVnJNwWoJLQqqEhiOomTfFMtvjQohba4Etbf04KERddxKcaaNWu6XQuoiv6SiFbUFZGC0UmcZEIr6iagxQVFgqdrASUfKiZA6YJCkcLTCaHkomICVKIl3OcvSLg6iIhkqfMJ4TkMGjKEa66dpRNCyUhFzEGpKoRIYcWrMLFs9nSqqw7jN8/9mg8//FBzVhVMc1AxVBVCpHASVZiYd/t9vPrqK4w/9TQNtUtaKiJA6YLC8hRbRDb2Vl9bW+ymVbT2E8K6uCeE1UfV8d4Hf9FQu6SlIgKULigsTwmXf49zzZQUzuDBg3k30WrTzXuZs2KlavdJWioiQEHyqhDKNBIJTnV1NZdfdjnLZk/vdEJ454Jr+dteNQwafkyn7TXULolUTICKXcL9iXWPEtm5g1uXLmH+goUaDxcJ2HdWLKe66jC+Pul0Zp59GrO/9Hn6DRnKB++3aqhd0lYxAapD7AWFSj0XyY+qqip+89yvuWraVbz/bjNHHdWb9U88yugxYzTULmmriDTzeJR6Xvq0kGFpiC2D1LNnT9Xuk7TTzCv2N0Jr2ZS+upoaLEERWQmPjlGLDiuWL2PxzYt0HZSkVLEBSgsSlj5VMi9dXYOWSDwVNwfVQannIiLhFkiAMrNzzGyzmW0zs24ZBmZ2sZm9Er39j5mNCeK4udKChCIi4ZVzkoSZHQ5sAc4GIsALwIXu/nrMNqcDG9292czOBRa5+ymp9l2o9aASrWWjNW5ERIJXyFp844Ft7r7d3Q8ADwJTYzdw9/9x9+bo3d8CDQEcNzBd17LR0hwiIsUXRIAaAuyMuR+JPpbINODRRE+a2XQzW29m6/fs2ZNz47KpEqHro0REii+IAGVxHos7bmhmn6c9QC1ItDN3v8/dG929sV+/flk3KtteUKJKzKoXJiJSWEEEqAgwNOZ+A9CtqJaZjQZWAVPdfW8Ax00q216QluYQEQmHIALUC8CxZnaMmR0BXACsjd3AzIYBPwEucfctARwzqVx6QVqaQ0QkHHIOUO7eBlwNrAM2Ag+5+wYzm2FmM6Kb3Qj0Ae4xsyYzy2tqXi69IF0fJSISDoFUknD3R4BHujy2Mubnq4CrgjhWOnKtEtFR5XzulInd6oVJ+NTX1sZdA6oHcDDO9nU1NapCIVICyrZY7Jy583jm+fWHhvk6ekFnjm9kxfJlae1D10GVhqRFYxM9HuLfe5Fyl+51UGUboNra2lQ1uUIoQImUlooPUB3UCyp/ClAipUXLbUSparKISGmq2GrmIiISbgpQUvLqamow6HbrEecxQwsahkE2Jcik8ihAScnbt38/7t7tdiDOY+6uFPMiUiFmyYQClIgUjAoxSyYUoESkIFpbW1m1ahX/9LXr6NmrPaNWhZglGQUoEcm7trY2rp01i7/85S/cdcN1fPXzjdx/62I+amtTIWZJSAFKRPJu/oKFNG3exj2//A13r3uOO37+NG9uep0frLhFhZglobK/DkpEiqtjdYEVP/tVp9UFrllyB7P+YQJvbnxNhZglLvWgRCSvkq0u0OOIv+GkTx2bViFmpaZXHgUoEcmrZGusfdx2gO/eeWfS+phKTa9cClAiklfJ1libduW0lEN7Sk2vXGVfLFZEii/b1QVaW1sZ0jC00/wVtAe4uVMmEtm5Q3NXJUjFYtOgSucihVFVVcWK5ctYfPOijP7m0lkdW8Wgy1dFDvFpTFukODpWF0j3hDDZ/JVS08tfRQYojWmLlIZk81dKTU9PfW0tZtbtVl9bm5fXBSmQOSgzOwe4AzgcWOXuS7o8b9HnJwN/Bi5399+n2m8+5qCSjWlf96UJbNm0iYEDBwZ6TAmX+tpamltauj1eV1OjQrIhpNWxc5N0Qc8k3//Zvi7NNqU1B5VzD8rMDgfuBs4FRgEXmtmoLpudCxwbvU0H7s31uNlKNqZddURPPnnsSA33lZmuZ4LNLS04dLvFC1pSfB3zV5GdO3hi3aNEdu5gxfJlCk4VIIghvvHANnff7u4HgAeBqV22mQp839v9FuhtZoMCOHbGko1pHzzwF5Y+/JiG+8pM14AkpSnT+SspfUEEqCHAzpj7kehjmW4DgJlNN7P1ZrZ+z549ATSvs0Rj2t9dOIuJ/3wBg4Yfo+rKIiIhEESAsjiPdT1RTWeb9gfd73P3Rndv7NevX86Ni+fWpUs4c3wj131pAld99iRmf+nzHP2pUXxlzg0Aqq4sIhICQQSoCDA05n4D0PWbPZ1tCqZjTHvLpk0c+ODPLPnRL7hs/o0cHh3TVgqriJSLupoaDLrdekDSLL1Er6urqSlY24MIUC8Ax5rZMWZ2BHABsLbLNmuBS63dqcB77v5WAMfOycCBA5k+fTprvnWDUlgrSB3d/+gK/YcnUij79u/H3bvdDtI9USg2WSjR6wqZ6ZpzGoy7t5nZ1cA62tPM17j7BjObEX1+JfAI7Snm22hPM78i1+MG5dalS5i/YCFzp0zslsIq5aGupgZTWrlIyVEtviiVPRKRSpLP65zSOLZq8WWiI4VVRETCoSJLHYmISPgpQEnZCUMNMZGwC0OWXioa4pOy01E5oqt4iRIilaoUEoTUgxIRkVBSgBIRkVBSgBJJQXNaIsWhOSiRFDSnJVIc6kFJ2SmF7CQRSU09KCk7pZCdJCKpqQclIiKhpAAlIiKhpAAlkoLmtESKQwFKJIUwrIsjko5yuyRCAUokgXL7Y5fy13FJRKJFCEuNsvhEEtD1TyLFpR6UiIiEkgKUSA40DCiSPzkFKDOrN7NfmtnW6L91cbYZama/MrONZrbBzGblckyRMCm3MX+pPGE+ycq1B7UQeNLdjwWejN7vqg2Y6+7HA6cCXzezUTkeV0REusjmkogwn2TlGqCmAvdHf74f+HLXDdz9LXf/ffTnFmAjMCTH44rkna5/EimuXLP4Brj7W9AeiMysf7KNzWw4cBLwuxyPK5J3us5JSk25ZZ6mDFBm9gQwMM5T38zkQGZWDfwYmO3uCf/yzWw6MB1g2LBhmRxCRETKSMoA5e5fSPScmb1tZoOivadBwO4E2/WgPTj9p7v/JMXx7gPuA2hsbIx3MiASGlW0D/vFe1wk3+pra7vNFRlQB+wrSouClesc1FrgsujPlwE/7bqBmRmwGtjo7ityPJ5IqNQkmI9K9LhIkBImOGSwjzDPteYaoJYAZ5vZVuDs6H3MbLCZPRLd5gzgEmCimTVFb5NzPK5IKCSq0weENnVXJFaYa03mNBLh7nuBs+I8vguYHP3518QfBREpW+U2WS2lJ/ZLNwy9oWxoqFxEpAx19ORLmUodiYhIKClAiYiUqDAnOARBQ3wiIiUqDIkM+aQelEgelPuZrWQnzIVZw0g9KJE8KPczW8lOS4IszkSPVzr1oERCQmfX5a+N7hfVevRx6U49KJEiiVemRtdOifyVApRIkXS9mFdXs4t0piE+EREJJQUoEREJJQUoEZEC0eUHmdEclEhI1BF/HkoWsL4pAAAIC0lEQVRfXuVDlx9kRgFKpEjqamriZujV1dToi0wEBSiRolEQEklOc1AiIkWii7OTU4ASEUlDPoJJwiXbdXE2oCE+EZG0aJXkwlMPSkREQimnAGVm9Wb2SzPbGv23Lsm2h5vZS2b281yOKVLpsh1q0nyHlJpce1ALgSfd/Vjgyej9RGYBG3M8nkjFy3beQvMdUmpyDVBTgfujP98PfDneRmbWAPw9sCrH44mIlA1Vlkgu1ySJAe7+FoC7v2Vm/RNsdzswH0j5v25m04HpAMOGDcuxeSIiwUh2YXW2dC1ccil7UGb2hJm9Fuc2NZ0DmNk/ALvd/cV0tnf3+9y90d0b+/Xrl85LRMqG5okyV6j/s3379+Pu3W6Jgow+y9yl7EG5+xcSPWdmb5vZoGjvaRCwO85mZwBTzGwy0BOoNbMfuPtXsm61SJlSKnPmwvp/FtZ2lZJc56DWApdFf74M+GnXDdz9endvcPfhwAXAUwpOItnLdt5C8x1SanKdg1oCPGRm04AdwPkAZjYYWOXuk3Pcv4h0ke28heY7pNTkFKDcfS9wVpzHdwHdgpO7Pw08ncsxRUSkMqiShIiIhJIClEiIFHueqBQzz4r9f1Zq7SolKhYrEiLZzhPV19bGrQiR6eKH8TLP6qOPm3Ve7zcsCyuGoQ3xhLVdpUQ9KJEykM8yRs1x9hvGEkml2PuT5NSDEpGyoOuOyo96UCIiEkoKUCJSFBqSk1QUoETkkHiZZ/lSyOU/FPRKkwKUSBkIKqU5XkHUINOlY3tNxOyrPuM9ZaZr0FPvrTQoSUKkDOQzpTnIfSdMZAhg3wmXw6A9EzGtdiihIlTUgxKRQBWrd9IRSLsOGe5L1V46D2eqNxUe6kGJSKBKrXfScZ1XV2FtbyVRgBKRUOjowagUkHRQgBKRUHCP14/JTD6WZZfiUYASkYLJdwBJN6EjUTskXBSgRKRgwlJANbYdXYvgSngoi09EAlVqy0yUWnsriXpQIhKosPSS0lVq7a0k6kGJiEgo5RSgzKzezH5pZluj/9Yl2K63mT1sZpvMbKOZnZbLcUVEpPzl2oNaCDzp7scCT0bvx3MH8Ji7fwoYA2zM8bgiIlLmcg1QU4H7oz/fD3y56wZmVgt8DlgN4O4H3P3dHI8rIiJlLtcANcDd3wKI/ts/zjYjgD3A98zsJTNbZWa9Eu3QzKab2XozW79nz54cmyciIqUqZYAysyfM7LU4t6lpHqMKGAfc6+4nAe+TeCgQd7/P3RvdvbFfv35pHkJERMpNyjRzd/9CoufM7G0zG+Tub5nZIGB3nM0iQMTdfxe9/zBJApSIiAjkPsS3Frgs+vNlwE+7buDufwJ2mtlx0YfOAl7P8bgiIlLmcg1QS4CzzWwrcHb0PmY22MweidnuGuA/zewVYCxwS47HFRGRMmdBVBDOFzPbA7yZp933Bd7J077DQu+xPOg9lge9x7862t1TJhmEOkDlk5mtd/fGYrcjn/Qey4PeY3nQe8ycSh2JiEgoKUCJiEgoVXKAuq/YDSgAvcfyoPdYHvQeM1Sxc1AiIhJuldyDEhGREFOAEhGRUKqYAGVm55vZBjP72MwSpkGa2TlmttnMtplZSZVkymB9rjfM7FUzazKz9YVuZzZSfS7W7s7o86+Y2bhitDMXabzHCWb2XvRzazKzG4vRzmyZ2Roz221mryV4vhw+w1TvsaQ/QwAzG2pmv4qu7bfBzGbF2SaYz9LdK+IGHA8cBzwNNCbY5nDgD7RXYD8CeBkYVey2Z/AebwUWRn9eCCxNsN0bQN9itzeD95XycwEmA48CBpwK/K7Y7c7De5wA/LzYbc3hPX6O9sLRryV4vqQ/wzTfY0l/htH3MAgYF/25BtiSr7/HiulBuftGd9+cYrPxwDZ33+7uB4AHaV/zqlSkXJ+rRKXzuUwFvu/tfgv0jhYwLhWl/ruXkrs/C+xLskmpf4bpvMeS5+5vufvvoz+30L4A7ZAumwXyWVZMgErTEGBnzP0I3f/jwyyd9bkAHHjczF40s+kFa1320vlcSv2zS7f9p5nZy2b2qJmdUJimFUypf4bpKpvP0MyGAycBv+vyVCCfZcrlNkqJmT0BDIzz1DfdvVul9Xi7iPNYqPLwk73HDHZzhrvvMrP+wC/NbFP0zC+s0vlcQv/ZpZBO+39Pew2zVjObDPw/4Ni8t6xwSv0zTEfZfIZmVg38GJjt7vu7Ph3nJRl/lmUVoDzJ2lVpigBDY+43ALty3Gegkr3HNNfnwt13Rf/dbWb/TfvwUpgDVDqfS+g/uxRStj/2S8DdHzGze8ysr7uXSwHSUv8MUyqXz9DMetAenP7T3X8SZ5NAPksN8XX2AnCsmR1jZkcAF9C+5lWpSLk+l5n1MrOajp+BLwJxM45CJJ3PZS1waTR76FTgvY7hzhKR8j2a2UAzs+jP42n/+91b8JbmT6l/himVw2cYbf9qYKO7r0iwWSCfZVn1oJIxs38Evgv0A35hZk3uPsnMBgOr3H2yu7eZ2dXAOtqzqta4+4YiNjtTS4CHzGwasAM4H9rX5yL6HoEBwH9H/0aqgP9y98eK1N60JPpczGxG9PmVwCO0Zw5tA/4MXFGs9mYjzfd4HvA1M2sDPgAu8GjKVCkwsx/SnsXW18wiwE1ADyiPzxDSeo8l/RlGnQFcArxqZk3Rx24AhkGwn6VKHYmISChpiE9EREJJAUpEREJJAUpEREJJAUpEREJJAUpEREJJAUpEREJJAUpERELp/wMD6aH7VhLRlQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[y_db == 0, 0], X[y_db == 0, 1],\n",
    "            c='lightblue', marker='o', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 1')\n",
    "plt.scatter(X[y_db == 1, 0], X[y_db == 1, 1],\n",
    "            c='red', marker='s', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 2')\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Accelerating DBSCAN with RAPIDS\n",
    "\n",
    "While the world’s data doubles each year, CPU computing has hit a brick wall with the end of Moore’s law. For the same reasons, scientific computing and deep learning has turned to NVIDIA GPU acceleration, data analytics and machine learning where GPU acceleration is ideal. \n",
    "\n",
    "NVIDIA created RAPIDS – an open-source data analytics and machine learning acceleration platform that leverages GPUs to accelerate computations. RAPIDS is based on Python, has pandas-like and Scikit-Learn-like interfaces, is built on Apache Arrow in-memory data format, and can scale from 1 to multi-GPU to multi-nodes. RAPIDS integrates easily into the world’s most popular data science Python-based workflows. RAPIDS accelerates data science end-to-end – from data prep, to machine learning, to deep learning. And through Arrow, Spark users can easily move data into the RAPIDS platform for acceleration.\n",
    "\n",
    "So how do we use RAPIDS? First, we cast our data to a `pandas.DataFrame` and use that to create a `cudf.DataFrame`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import cudf\n",
    "\n",
    "X_df = pd.DataFrame({'fea%d'%i: X[:, i] for i in range(X.shape[1])})\n",
    "X_gpu = cudf.DataFrame.from_pandas(X_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we load the `DBSCAN` class from the `cuml` package and instantiate it in the same way we did with the `sklearn.cluster.DBSCAN` class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cuml import DBSCAN as cumlDBSCAN\n",
    "\n",
    "db_gpu = cumlDBSCAN(eps=0.2, min_samples=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `DBSCAN` class from `cuml` implements the same API as the `sklearn` version; we can use the `fit` and `fit_predict` methods to fit our `DBSCAN` model to the data and generate predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_db_gpu = db_gpu.fit_predict(X_gpu)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lastly, let's visualize our results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XuUVOWV9/Hv1sbBobulmzs0iCRiRAeQtHhLIsEYlJlAZkbf5SVecQghKggE0KylyCQuUCBqvDAOkDFmJsbXZN6QRMWoUVecJIqxvSDXMAoVjCC00m000LrfP7qaVHfXvU5Vnar6fdaqRVfVqXOeorprn+d59tmPuTsiIiJhc1ixGyAiIhKPApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiISSApSIiIRSVbEbkEzfvn19+PDhxW6GiIgE6MUXX3zH3ful2i7UAWr48OGsX7++2M0QEZEAmdmb6WynIT4REQklBSgREQklBSgREQmlUM9BiYgU28GDB4lEInz44YfFbkrJ6dmzJw0NDfTo0SOr1ytAiYgkEYlEqKmpYfjw4ZhZsZtTMtydvXv3EolEOOaYY7Lah4b4RESS+PDDD+nTp4+CU4bMjD59+uTU81SAEhFJQcEpO7n+vylASUVpbW1ly5YttLa2FrspIpKCApSEXhBBpa2tjTlz5zGkYShnTTqXIQ1DmTN3Hm1tbQG2VKRdIU6EFi1axLJlyzJ+3bvvvss999yT8/HvuusuPvnJT2JmvPPOOznvLx4FKAmtIIPK/AULeeb59az42a+447Ffs+Jnv+KZ59czf8HCPLRcKlUpnAhlE6DcnY8//rjTY2eccQZPPPEERx99dJDN60QBSkIrqKDS2trK6tWrmXnL7dT1HwBAXf8BzLzldtasWdPpLFdDgJKLfJ0Iff/732f06NGMGTOGSy65pNvzEyZMOFQW7p133qGjhumGDRsYP348Y8eOZfTo0WzdupWFCxfyhz/8gbFjx/KNb3wDgNtuu42TTz6Z0aNHc9NNNwHwxhtvcPzxxzNz5kzGjRvHzp07Ox3zpJNOIt+1UhWgJJQyCSqp7Nq1i9q6+kP76VDXfwC1vevYtWtXSZz5SrgF+Tsba8OGDXz729/mqaee4uWXX+aOO+5I+7UrV65k1qxZNDU1sX79ehoaGliyZAmf+MQnaGpq4rbbbuPxxx9n69atPP/88zQ1NfHiiy/y7LPPArB582YuvfRSXnrppbz2lBJRgJJQSieopGvw4MHsb95H8+63Oz3evPtt9r/bzODBgzUEKDkL8nc21lNPPcV5551H3759Aaivr0/7taeddhq33HILS5cu5c033+TII4/sts3jjz/O448/zkknncS4cePYtGkTW7duBeDoo4/m1FNPzardQVCAklBKJ6ikq7q6mmnTpnHPDbMP7a9599vcc8NsrrzySoC8nPlKZQnydzaWu6dM166qqjo0RxR73dFFF13E2rVrOfLII5k0aRJPPfVU3P1ff/31NDU10dTUxLZt25g2bRoAvXr1yqrNQVGAklBKFVSqq6sz2t+tS5dw5vhG5k6ZyOxzPsPcKRM5c3wjty5dkrcz3yBoTqx0BP072+Gss87ioYceYu/evQDs27ev2zbDhw/nxRdfBODhhx8+9Pj27dsZMWIE1157LVOmTOGVV16hpqaGlpaWQ9tMmjSp04nYH//4R3bv3p1VW4MWSIAyszVmttvMXkvwvJnZnWa2zcxeMbNxQRxXSk+6X7itra1cNe1KTj9pTNygkqmqqipWLF9GZOcOnlj3KJGdO1ixfBlVVVV5O/PNhebESlOyE6FsnXDCCXzzm9/kzDPPZMyYMcyZM6fbNvPmzePee+/l9NNP75Ty/aMf/YgTTzyRsWPHsmnTJi699FL69OnDGWecwYknnsg3vvENvvjFL3LRRRdx2mmn8Xd/93ecd955nQJYInfeeScNDQ1EIhFGjx7NVVddlfV7TMjdc74BnwPGAa8leH4y8ChgwKnA79LZ76c//WmX8nDw4EG/bs5crz2qtzcMH+G1R/X26+bM9YMHD6bc7utXX+MbNmzwlpaWvLXvujlzfdxnzvRVz77kP960y1c9+5KP+8yZft2cuSlf29LS4ps3bw60fbm0R4L1+uuvZ/yafPxOlKp4/3/Aek8ntqSzUVo7guFJAtS/ARfG3N8MDEq1TwWo8pHuF26xvpg7AuNRvet86PARflTvurgBNN5rUgXdTLW0tHjtUb0P/R903FY9+5If1btOX3oFlk2Akr/KJUAVag5qCBCbRB+JPtaNmU03s/Vmtn7Pnj0FaZzEF9T8R7rpt/lI0033PSQbAkwkX5l/YZ4TEymkQgWoeCkoHm9Dd7/P3RvdvbFfv355bpbEE/T8R7pfuEF+MWf7Hqqrqxk5cmTKCe18XfMC+csGEyk1hQpQEWBozP0GQKeBIRV0zyDdL9wgv5jzfV1TPns5+coGEyk1hQpQa4FLo9l8pwLvuftbBTq2pKFjKOxPf/pT4D2DdL9wg/pizmfvpkO+ezn5yAYTKTVBpZn/EPgNcJyZRcxsmpnNMLMZ0U0eAbYD24B/B2YGcVzJXdehsGOPO47DevSgtr5Pp+1y7Rmk+4UbxBdzPno3Xeey8t3LyWZOTKTspJNJUaybsvjyL17W3KfGneyTLrg0Lxlk6abf5pKmG2QWXLJMvWwy/6T0hDWL76abbvLbbrst49c1Nzf73XffnfPxL7roIh85cqSfcMIJfsUVV/iBAwfiblcKWXwSQomGwubdfh9P//T/8tYb/wsE2zNINwkh3e0SvTao3k2yuax4vZzFNy9i+/btnTITVQmictTX1mJm3W71tbXFbtohQS23cfHFF7Np0yZeffVVPvjgA1atWhVkMwGVOqpoO3bs4G/+tleCobB6Fp5/bsnOfwQxVJjuXFZ1dTUjRozgxpsWxWQNNnDyKacyeEiDKkFUkOaWFhy63ZrTqMyQTBiX25g8efKhADx+/HgikUhO7zEeDWhXsHvuXXlooj82SDXvfpsDH/yZrVs2s3//fgYPHlxymWMdvZvFNy9i165dWb2HdOayRo4cCXTuadX1H0Dz7rdZPmcGp0+eyvRFSw714OYvWMiK5ZmvgpqN1tbWrN+7hEfHchvPPfccffv2jVuLL5GO5TYuvvhiDhw4wEcffcSSJUt47bXXaGpqAui03Ia7M2XKFJ599lmGDRvG5s2b+d73vpe0x3Xw4EEeeOCBjJYBSZd6UBWqtbWVBx54gAlfPp/vLpzVaShs2ezpXHLJJQwcODDrYbawyGWoMN1MvUQ9rbkrVvLco2v54P33C1odXXX8ykvYl9uYOXMmn/vc5/jsZz+bwbtKjwJUheroHVx03fX0HTyEWf8wga9POoPZX/o8b+94g6/N+Gqxm1h06c5lJetpVfeuo3n3nw7dL0QlCK1tVV7cw7vcxs0338yePXtYsWJFpm8rLQpQFap///7se2cPX/vCKWx44be4O6MaT+Hb//VTPm5rY9iwYcVuYiikM5eVrKfV+m4zdf0HHrqf70oQhbgGTAorrMttrFq1inXr1vHDH/6Qww7LTyhRgKpQi//1W4wYdSJ3/uIZ7l73HHf+4hne2fVHls68XNUKYqRzPVKintbyOTM449wpHNmrV8EqQaiOX3HV1dRg0O1WV1OT9T7DutzGjBkzePvttznttNMYO3Ysixcvzvo9JmLtKenh1NjY6B2ZKRKc1tZWhjQMPTSh36F599tcc+5n2BWJ0Lt37yK2sPS0tbUxf8FC1qxZQ23vOt57dx8jRx7Hls2bOaqunv3vNnPllVdy69Ileb3YNtlnO3fKRCI7d+jkI0MbN27k+OOPL3YzSla8/z8ze9HdG1O9Vll8FSjZWXbffgPYvXu3AlSGEmUNFjqTLrY31zHMpzp+UqoUoCpQ7JxJ17PsQlfLLrdU6I6swUT3C+HWpUuYv2Ahc6dMpLZ3Xafem0gp0RxUBQpDtWylQueP6vgFL8xTIWGW6/+bAlSFKna1bKVC518u14DJX/Xs2ZO9e/cqSGXI3dm7dy89e/bMeh9KkqgAyYbRijHEpol8KSUHDx4kEol0ur5I0tOzZ08aGhro0aNHp8eVJCGHMstWr15NbV09+5v3MW3atE6ZZMWYI8mkhJBIsfXo0YNjjjmm2M2oSBriK2NhHUbTkuYikg4FqDIV5ooCYUjSkOS0TIiEgQJUmQp7RYFiJ2lIfMqulDBRgCpTYR9GUyp0OIV1WFgqUyAByszOMbPNZrbNzLr9JpvZUWb2MzN72cw2mNkVQRy30iUbhimVYTSlQodHmIeFpTLlHKDM7HDgbuBcYBRwoZmN6rLZ14HX3X0MMAFYbmZH5HrsSpXuMIyG0SQTYR8WlsoTRA9qPLDN3be7+wHgQWBql20cqLH2RU2qgX2ABrWzlO4wjIbRJBNhHxaWyhNEgBoCxC5WH4k+Fusu4HhgF/AqMMvdPw7g2BUnm2EYDaNJOkplWFgqRxABKt5Sj13LU0wCmoDBwFjgLjOrjbszs+lmtt7M1u/ZsyeA5pUXDcNUrkKkfmtYWMIkiAAVAYbG3G+gvacU6wrgJ95uG/C/wKfi7czd73P3Rndv7NevXwDNKy8ahqk8hUz91rCwhEkQAeoF4FgzOyaa+HABsLbLNjuAswDMbABwHLA9gGOXvEzPijUMU3mKkfqtYWEJg0CKxZrZZOB24HBgjbt/28xmALj7SjMbDPwHMIj2IcEl7v6DVPst52Kx6dTJS/XajtVbC7VaqxReroV1y229LSkP6RaLVTXzIpkzdx7PPL++26qnZ45vZMXyZWntQ18+5W/Lli2cNelc7njs192em33OZ3hi3aNxC+vmcgIkkm/pBihVkiiCoC6I1DBM+ct2zlEVIaQcKEAVQUcmXs9e1ez63z/wwfvvA8rEk+6ymXNURQgpF+rrF0H//v3Z984epk/4NLX1fWhp3sdZ513I339lmjLxpJtbly5h/oKFzJ0ysducYzxab0vKhQJUESz+128xYtSJzF52z6H5pzsXXMvNV/wfZeJJNx2p34tvXpTWnGPssGDXxAqdAEkpUZJEgSXLyrrm3M+wKxKhd+/eRWyhlINMknCUbCOFpiSJkEo2/NK33wB2796d0f60sJzEk05FCK39JGGnAFVgQVWC0JeLJJNORQhl+knYKUAVWFCVIPTlIulIdCmCMv2kFChAFUGuBTn15SK5UtFhKQUKUEWQa0FOfblIrlR0WEqBAlQRZVsJQl8ukisVHZZSoABVgvTlIkHQ2k8SdroOqkSporkERddBSaGpmnmF0JeLiJQaXahbIVTRXHKli70lrBSgRCpUqou9Fbik2DRZIVKhYi/2jq3XN2/+fA6zw7TYoRSd5qBEKlCqosUjR5/ENUu/m/VqzyLJaA5KRBJKdrH33xzZi/OvnqcqJVJ0gQQoMzvHzDab2TYzi1sMzswmmFmTmW0ws2eCOK6IZCfZxd5/bt3P8E+d0OlxVSmRYsg5QJnZ4cDdwLnAKOBCMxvVZZvewD3AFHc/ATg/1+OKSPYSXex99/WzqKrqwYfvd+4pqUqJFEMQPajxwDZ33+7uB4AHgaldtrkI+Im77wBw98wWPSphyoSSsIpXSWLCKSfzL1ddpSolEgpBpOQMAXbG3I8Ap3TZZiTQw8yeBmqAO9z9+/F2ZmbTgekAw4YNC6B5xdFR6UGZUBJWiZaS7/jdnTtlYrcqJSKFlHMWn5mdD0xy96ui9y8Bxrv7NTHb3AU0AmcBRwK/Af7e3bck23epZfHFVnW48aZFaS+5LRJGqlIi+ZJuFl8Qp/IRYGjM/Qag60xqBHjH3d8H3jezZ4ExQNIAVSq69pbea97LRx99zHfWPtUtE2rulIksvnmR/uAl9DqqlGRKgU2CEsQc1AvAsWZ2jJkdAVwArO2yzU+Bz5pZlZn9Le1DgBsDOHYodF3d9js/e5phxx3PL36wutN2yoSScpaqMoVIpnIOUO7eBlwNrKM96Dzk7hvMbIaZzYhusxF4DHgFeB5Y5e6v5XrsMEi0uu287/wbTz78Qz54//1D2yoTSspZ1xO1FT/7Fc88v575C+JeeSKSUiDXQbn7I+4+0t0/4e7fjj620t1Xxmxzm7uPcvcT3f32II4bBskveDySNzZtAJQJJeUt0YmaLvCVXCidLEexFzx2LRnzQWsLt868nKPq6pUJJWUt2Ylax7B2NvNZUtlU6ihHyVa3nfHVGfwxspMn1j1KZOcOVixfphRzKUvJKlNoWFuypQAVgGRLZ2u9JqkEyU7UNKwt2VI18wApvVYqWcflFmvWrOl2ga9GDiSWlnwXkaKIPVEDdNIm3Wi5DREpiurqakaMGMGNNy3SNVGSEwUoEQmcromSIGiIT0QClWy13rlTJhLZuUPDfRVOQ3wiUhTpXBMlkg4FKBEJlK6JkqAoQIlIRlItwqlroiQoClBSkuprazGzbrf62tpiN61sZVKtPNnF6yLpUpKElCQzI95vrgFh/p0uZXPmzst4EU5dvC7x6ELdgOgPLJxiA1Q90Bxnm7qaGvbt31/AVpUvZeZJkJTFlyMtvlY6mgGPc2tuaSlms8qKMvOkGBSgEkh0oeF1c+Z2myBONWksUuqUmSfFoCG+OJINZ8w8+1T6DhhE63vvcsUVV4DB99Z8j9q6evY372PatGkqjlkAsUN8BpqPKoBs5qBE4kl3iE/fonEkG86o6z+Q6//tAY7sVc2Nl/wT9QMGHgpkHX+w8xcs1B9sntXV1GAawiuoW5cuYf6ChcydMrFbtXKRfFAPKo5kPajZX/o8K596AYDpEz7Nnb94RpPGBVRfWxt3bqkO2NflMfWg8kOJQ5KrgiZJmNk5ZrbZzLaZWcJqkGZ2spl9ZGbnBXHcfEl0oeF3F85i4j9fwJG9etG8+0/U1vfRpHGBNbe0xE+IoD0gxd7qamqK1cyypkU4pVByDlBmdjhwN3AuMAq40MxGJdhuKbAu12MWQuyFhrPOOYOZZ5/KwGHD+cqcGwCo6z+Q/fv2atI4RNy9000p5iKlLYge1Hhgm7tvd/cDwIPA1DjbXQP8GNgdwDHzrqqqihXLlxHZuYMn1z3Gv/zLdPZE3mT/vr0AfPh+K73r+3D7vJkq5yIikgdBJEkMAXbG3I8Ap8RuYGZDgH8EJgInJ9uZmU0HpgMMGzYsgOblprq6msGDBzPzazPg3pWdJogvv/xyMDRpLCKSB0EEKIvzWNeZ6duBBe7+kVm8zWNe6H4fcB+0J0kE0L6stbW1MX/BQlavXn0ojfzSSy/lazO+yrBhww71kr61eLEmjaWiKXFC8iGIIb4IMDTmfgPQNUOgEXjQzN4AzgPuMbMvB3DsvIp3se7/vPQyq1av6fRHqEnj/IotDFtF92QIJUQUjyquSD7lnGZuZlXAFuAs4I/AC8BF7r4hwfb/Afzc3R9Ote9i1uJT7bHwUGHY8NLFu5KNgqWZu3sbcDXt2XkbgYfcfYOZzTCzGbnuv1hUe0wkudbWVlavXn0oOEH738fMW25nzZo1Kv0lOQukkoS7PwI80uWxlQm2vTyIY+ZbbO2xrj0opZGLpHcSN3LkyCK1TsqBisUmoFVBRZJTAVnJNwWoJLQqqEhiOomTfFMtvjQohba4Etbf04KERddxKcaaNWu6XQuoiv6SiFbUFZGC0UmcZEIr6iagxQVFgqdrASUfKiZA6YJCkcLTCaHkomICVKIl3OcvSLg6iIhkqfMJ4TkMGjKEa66dpRNCyUhFzEGpKoRIYcWrMLFs9nSqqw7jN8/9mg8//FBzVhVMc1AxVBVCpHASVZiYd/t9vPrqK4w/9TQNtUtaKiJA6YLC8hRbRDb2Vl9bW+ymVbT2E8K6uCeE1UfV8d4Hf9FQu6SlIgKULigsTwmXf49zzZQUzuDBg3k30WrTzXuZs2KlavdJWioiQEHyqhDKNBIJTnV1NZdfdjnLZk/vdEJ454Jr+dteNQwafkyn7TXULolUTICKXcL9iXWPEtm5g1uXLmH+goUaDxcJ2HdWLKe66jC+Pul0Zp59GrO/9Hn6DRnKB++3aqhd0lYxAapD7AWFSj0XyY+qqip+89yvuWraVbz/bjNHHdWb9U88yugxYzTULmmriDTzeJR6Xvq0kGFpiC2D1LNnT9Xuk7TTzCv2N0Jr2ZS+upoaLEERWQmPjlGLDiuWL2PxzYt0HZSkVLEBSgsSlj5VMi9dXYOWSDwVNwfVQannIiLhFkiAMrNzzGyzmW0zs24ZBmZ2sZm9Er39j5mNCeK4udKChCIi4ZVzkoSZHQ5sAc4GIsALwIXu/nrMNqcDG9292czOBRa5+ymp9l2o9aASrWWjNW5ERIJXyFp844Ft7r7d3Q8ADwJTYzdw9/9x9+bo3d8CDQEcNzBd17LR0hwiIsUXRIAaAuyMuR+JPpbINODRRE+a2XQzW29m6/fs2ZNz47KpEqHro0REii+IAGVxHos7bmhmn6c9QC1ItDN3v8/dG929sV+/flk3KtteUKJKzKoXJiJSWEEEqAgwNOZ+A9CtqJaZjQZWAVPdfW8Ax00q216QluYQEQmHIALUC8CxZnaMmR0BXACsjd3AzIYBPwEucfctARwzqVx6QVqaQ0QkHHIOUO7eBlwNrAM2Ag+5+wYzm2FmM6Kb3Qj0Ae4xsyYzy2tqXi69IF0fJSISDoFUknD3R4BHujy2Mubnq4CrgjhWOnKtEtFR5XzulInd6oVJ+NTX1sZdA6oHcDDO9nU1NapCIVICyrZY7Jy583jm+fWHhvk6ekFnjm9kxfJlae1D10GVhqRFYxM9HuLfe5Fyl+51UGUboNra2lQ1uUIoQImUlooPUB3UCyp/ClAipUXLbUSparKISGmq2GrmIiISbgpQUvLqamow6HbrEecxQwsahkE2Jcik8ihAScnbt38/7t7tdiDOY+6uFPMiUiFmyYQClIgUjAoxSyYUoESkIFpbW1m1ahX/9LXr6NmrPaNWhZglGQUoEcm7trY2rp01i7/85S/cdcN1fPXzjdx/62I+amtTIWZJSAFKRPJu/oKFNG3exj2//A13r3uOO37+NG9uep0frLhFhZglobK/DkpEiqtjdYEVP/tVp9UFrllyB7P+YQJvbnxNhZglLvWgRCSvkq0u0OOIv+GkTx2bViFmpaZXHgUoEcmrZGusfdx2gO/eeWfS+phKTa9cClAiklfJ1libduW0lEN7Sk2vXGVfLFZEii/b1QVaW1sZ0jC00/wVtAe4uVMmEtm5Q3NXJUjFYtOgSucihVFVVcWK5ctYfPOijP7m0lkdW8Wgy1dFDvFpTFukODpWF0j3hDDZ/JVS08tfRQYojWmLlIZk81dKTU9PfW0tZtbtVl9bm5fXBSmQOSgzOwe4AzgcWOXuS7o8b9HnJwN/Bi5399+n2m8+5qCSjWlf96UJbNm0iYEDBwZ6TAmX+tpamltauj1eV1OjQrIhpNWxc5N0Qc8k3//Zvi7NNqU1B5VzD8rMDgfuBs4FRgEXmtmoLpudCxwbvU0H7s31uNlKNqZddURPPnnsSA33lZmuZ4LNLS04dLvFC1pSfB3zV5GdO3hi3aNEdu5gxfJlCk4VIIghvvHANnff7u4HgAeBqV22mQp839v9FuhtZoMCOHbGko1pHzzwF5Y+/JiG+8pM14AkpSnT+SspfUEEqCHAzpj7kehjmW4DgJlNN7P1ZrZ+z549ATSvs0Rj2t9dOIuJ/3wBg4Yfo+rKIiIhEESAsjiPdT1RTWeb9gfd73P3Rndv7NevX86Ni+fWpUs4c3wj131pAld99iRmf+nzHP2pUXxlzg0Aqq4sIhICQQSoCDA05n4D0PWbPZ1tCqZjTHvLpk0c+ODPLPnRL7hs/o0cHh3TVgqriJSLupoaDLrdekDSLL1Er6urqSlY24MIUC8Ax5rZMWZ2BHABsLbLNmuBS63dqcB77v5WAMfOycCBA5k+fTprvnWDUlgrSB3d/+gK/YcnUij79u/H3bvdDtI9USg2WSjR6wqZ6ZpzGoy7t5nZ1cA62tPM17j7BjObEX1+JfAI7Snm22hPM78i1+MG5dalS5i/YCFzp0zslsIq5aGupgZTWrlIyVEtviiVPRKRSpLP65zSOLZq8WWiI4VVRETCoSJLHYmISPgpQEnZCUMNMZGwC0OWXioa4pOy01E5oqt4iRIilaoUEoTUgxIRkVBSgBIRkVBSgBJJQXNaIsWhOSiRFDSnJVIc6kFJ2SmF7CQRSU09KCk7pZCdJCKpqQclIiKhpAAlIiKhpAAlkoLmtESKQwFKJIUwrIsjko5yuyRCAUokgXL7Y5fy13FJRKJFCEuNsvhEEtD1TyLFpR6UiIiEkgKUSA40DCiSPzkFKDOrN7NfmtnW6L91cbYZama/MrONZrbBzGblckyRMCm3MX+pPGE+ycq1B7UQeNLdjwWejN7vqg2Y6+7HA6cCXzezUTkeV0REusjmkogwn2TlGqCmAvdHf74f+HLXDdz9LXf/ffTnFmAjMCTH44rkna5/EimuXLP4Brj7W9AeiMysf7KNzWw4cBLwuxyPK5J3us5JSk25ZZ6mDFBm9gQwMM5T38zkQGZWDfwYmO3uCf/yzWw6MB1g2LBhmRxCRETKSMoA5e5fSPScmb1tZoOivadBwO4E2/WgPTj9p7v/JMXx7gPuA2hsbIx3MiASGlW0D/vFe1wk3+pra7vNFRlQB+wrSouClesc1FrgsujPlwE/7bqBmRmwGtjo7ityPJ5IqNQkmI9K9LhIkBImOGSwjzDPteYaoJYAZ5vZVuDs6H3MbLCZPRLd5gzgEmCimTVFb5NzPK5IKCSq0weENnVXJFaYa03mNBLh7nuBs+I8vguYHP3518QfBREpW+U2WS2lJ/ZLNwy9oWxoqFxEpAx19ORLmUodiYhIKClAiYiUqDAnOARBQ3wiIiUqDIkM+aQelEgelPuZrWQnzIVZw0g9KJE8KPczW8lOS4IszkSPVzr1oERCQmfX5a+N7hfVevRx6U49KJEiiVemRtdOifyVApRIkXS9mFdXs4t0piE+EREJJQUoEREJJQUoEZEC0eUHmdEclEhI1BF/HkoWsL4pAAAIC0lEQVRfXuVDlx9kRgFKpEjqamriZujV1dToi0wEBSiRolEQEklOc1AiIkWii7OTU4ASEUlDPoJJwiXbdXE2oCE+EZG0aJXkwlMPSkREQimnAGVm9Wb2SzPbGv23Lsm2h5vZS2b281yOKVLpsh1q0nyHlJpce1ALgSfd/Vjgyej9RGYBG3M8nkjFy3beQvMdUmpyDVBTgfujP98PfDneRmbWAPw9sCrH44mIlA1Vlkgu1ySJAe7+FoC7v2Vm/RNsdzswH0j5v25m04HpAMOGDcuxeSIiwUh2YXW2dC1ccil7UGb2hJm9Fuc2NZ0DmNk/ALvd/cV0tnf3+9y90d0b+/Xrl85LRMqG5okyV6j/s3379+Pu3W6Jgow+y9yl7EG5+xcSPWdmb5vZoGjvaRCwO85mZwBTzGwy0BOoNbMfuPtXsm61SJlSKnPmwvp/FtZ2lZJc56DWApdFf74M+GnXDdz9endvcPfhwAXAUwpOItnLdt5C8x1SanKdg1oCPGRm04AdwPkAZjYYWOXuk3Pcv4h0ke28heY7pNTkFKDcfS9wVpzHdwHdgpO7Pw08ncsxRUSkMqiShIiIhJIClEiIFHueqBQzz4r9f1Zq7SolKhYrEiLZzhPV19bGrQiR6eKH8TLP6qOPm3Ve7zcsCyuGoQ3xhLVdpUQ9KJEykM8yRs1x9hvGEkml2PuT5NSDEpGyoOuOyo96UCIiEkoKUCJSFBqSk1QUoETkkHiZZ/lSyOU/FPRKkwKUSBkIKqU5XkHUINOlY3tNxOyrPuM9ZaZr0FPvrTQoSUKkDOQzpTnIfSdMZAhg3wmXw6A9EzGtdiihIlTUgxKRQBWrd9IRSLsOGe5L1V46D2eqNxUe6kGJSKBKrXfScZ1XV2FtbyVRgBKRUOjowagUkHRQgBKRUHCP14/JTD6WZZfiUYASkYLJdwBJN6EjUTskXBSgRKRgwlJANbYdXYvgSngoi09EAlVqy0yUWnsriXpQIhKosPSS0lVq7a0k6kGJiEgo5RSgzKzezH5pZluj/9Yl2K63mT1sZpvMbKOZnZbLcUVEpPzl2oNaCDzp7scCT0bvx3MH8Ji7fwoYA2zM8bgiIlLmcg1QU4H7oz/fD3y56wZmVgt8DlgN4O4H3P3dHI8rIiJlLtcANcDd3wKI/ts/zjYjgD3A98zsJTNbZWa9Eu3QzKab2XozW79nz54cmyciIqUqZYAysyfM7LU4t6lpHqMKGAfc6+4nAe+TeCgQd7/P3RvdvbFfv35pHkJERMpNyjRzd/9CoufM7G0zG+Tub5nZIGB3nM0iQMTdfxe9/zBJApSIiAjkPsS3Frgs+vNlwE+7buDufwJ2mtlx0YfOAl7P8bgiIlLmcg1QS4CzzWwrcHb0PmY22MweidnuGuA/zewVYCxwS47HFRGRMmdBVBDOFzPbA7yZp933Bd7J077DQu+xPOg9lge9x7862t1TJhmEOkDlk5mtd/fGYrcjn/Qey4PeY3nQe8ycSh2JiEgoKUCJiEgoVXKAuq/YDSgAvcfyoPdYHvQeM1Sxc1AiIhJuldyDEhGREFOAEhGRUKqYAGVm55vZBjP72MwSpkGa2TlmttnMtplZSZVkymB9rjfM7FUzazKz9YVuZzZSfS7W7s7o86+Y2bhitDMXabzHCWb2XvRzazKzG4vRzmyZ2Roz221mryV4vhw+w1TvsaQ/QwAzG2pmv4qu7bfBzGbF2SaYz9LdK+IGHA8cBzwNNCbY5nDgD7RXYD8CeBkYVey2Z/AebwUWRn9eCCxNsN0bQN9itzeD95XycwEmA48CBpwK/K7Y7c7De5wA/LzYbc3hPX6O9sLRryV4vqQ/wzTfY0l/htH3MAgYF/25BtiSr7/HiulBuftGd9+cYrPxwDZ33+7uB4AHaV/zqlSkXJ+rRKXzuUwFvu/tfgv0jhYwLhWl/ruXkrs/C+xLskmpf4bpvMeS5+5vufvvoz+30L4A7ZAumwXyWVZMgErTEGBnzP0I3f/jwyyd9bkAHHjczF40s+kFa1320vlcSv2zS7f9p5nZy2b2qJmdUJimFUypf4bpKpvP0MyGAycBv+vyVCCfZcrlNkqJmT0BDIzz1DfdvVul9Xi7iPNYqPLwk73HDHZzhrvvMrP+wC/NbFP0zC+s0vlcQv/ZpZBO+39Pew2zVjObDPw/4Ni8t6xwSv0zTEfZfIZmVg38GJjt7vu7Ph3nJRl/lmUVoDzJ2lVpigBDY+43ALty3Gegkr3HNNfnwt13Rf/dbWb/TfvwUpgDVDqfS+g/uxRStj/2S8DdHzGze8ysr7uXSwHSUv8MUyqXz9DMetAenP7T3X8SZ5NAPksN8XX2AnCsmR1jZkcAF9C+5lWpSLk+l5n1MrOajp+BLwJxM45CJJ3PZS1waTR76FTgvY7hzhKR8j2a2UAzs+jP42n/+91b8JbmT6l/himVw2cYbf9qYKO7r0iwWSCfZVn1oJIxs38Evgv0A35hZk3uPsnMBgOr3H2yu7eZ2dXAOtqzqta4+4YiNjtTS4CHzGwasAM4H9rX5yL6HoEBwH9H/0aqgP9y98eK1N60JPpczGxG9PmVwCO0Zw5tA/4MXFGs9mYjzfd4HvA1M2sDPgAu8GjKVCkwsx/SnsXW18wiwE1ADyiPzxDSeo8l/RlGnQFcArxqZk3Rx24AhkGwn6VKHYmISChpiE9EREJJAUpEREJJAUpEREJJAUpEREJJAUpEREJJAUpEREJJAUpERELp/wMD6aH7VhLRlQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(X[y_db_gpu == 0, 0], X[y_db_gpu == 0, 1],\n",
    "            c='lightblue', marker='o', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 1')\n",
    "plt.scatter(X[y_db_gpu == 1, 0], X[y_db_gpu == 1, 1],\n",
    "            c='red', marker='s', s=40,\n",
    "            edgecolor='black', \n",
    "            label='cluster 2')\n",
    "plt.legend()\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Benchmarking: Comparing GPU and CPU\n",
    "\n",
    "RAPIDS uses GPUs to parallelize operations and accelerate computations. We saw porting an example from the traditional scikit-learn interface to cuML was trivial. So how much speedup do we get from using RAPIDS? \n",
    "\n",
    "The answer to this question varies depending on the size and shape of the data. In the below example, we generate a matrix of 10,000 rows by 128 columns and show we were able to reduce computational time from ~45 seconds to ~5 seconds - almost a 9x speedup. Feel free to change the number of rows and columns to see how this speedup might change depending on the size and shape of the data.\n",
    "\n",
    "As a good rule of thumb, larger datasets will benefit from RAPIDS. There is overhead associated with using a GPU; data has to be transferred from the CPU to the GPU, computations have to take place on the GPU, and the results need to be transferred back from the GPU to the CPU. However, the transactional overhead of moving data back and forth from the CPU to the GPU can quickly become negligible due to the performance speedup from computing on a GPU instead of a CPU."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(10000, 128)\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "n_rows, n_cols = 10000, 128\n",
    "X = np.random.rand(n_rows, n_cols)\n",
    "print(X.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### GPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_df = pd.DataFrame({'fea%d'%i: X[:, i] for i in range(X.shape[1])})\n",
    "X_gpu = cudf.DataFrame.from_pandas(X_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "db_gpu = cumlDBSCAN(eps=3, min_samples=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 181 ms, sys: 39.1 ms, total: 220 ms\n",
      "Wall time: 219 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "y_db_gpu = db_gpu.fit_predict(X_gpu)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### CPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "db = DBSCAN(eps=3, min_samples=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 26.1 s, sys: 5.04 ms, total: 26.1 s\n",
      "Wall time: 26.1 s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "y_db = db.fit_predict(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Conclusion\n",
    "\n",
    "In conclusion, there are certain cases the DBSCAN algorithm can do a better job of clustering than traditional algorithms such as K Means or Agglomerative Clustering. Additionally, porting DBSCAN from CPU to GPU using RAPIDS is a trivial exercise and can yield massive performance gains.\n",
    "\n",
    "To learn more about RAPIDS, be sure to check out: \n",
    "\n",
    "* [Open Source Website](http://rapids.ai)\n",
    "* [GitHub](https://github.com/rapidsai/)\n",
    "* [Press Release](https://nvidianews.nvidia.com/news/nvidia-introduces-rapids-open-source-gpu-acceleration-platform-for-large-scale-data-analytics-and-machine-learning)\n",
    "* [NVIDIA Blog](https://blogs.nvidia.com/blog/2018/10/10/rapids-data-science-open-source-community/)\n",
    "* [Developer Blog](https://devblogs.nvidia.com/gpu-accelerated-analytics-rapids/)\n",
    "* [NVIDIA Data Science Webpage](https://www.nvidia.com/en-us/deep-learning-ai/solutions/data-science/)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Disclaimer: The above examples borrowed code snippets from *Python Machine Learning, 2nd Ed.* by Sebastian Raschka and Vahid Mirjalili. For a great deep dive into these concepts, the curious reader is strongly encouraged to explore that fantastic resource."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
